Bulls Eye for Athletes

Statistical analyses

Author
Affiliation

Gustaf Reinebo

Karolinska Institutet

Published

2025-02-06

1 Data preparation

The following data preparations are the basis for all analyses related to the Bulls Eye for Athletes (BEA) study. Code for scale step alterations and some of the modified data frames are presented next to their related sub-analysis to make the workflow more comprehensible.

Code
Code
#loading bulls eye for athletes dataset
BullsEye_merged_data_final_without_comp_level <- read_excel("BullsEye_merged_data_final.xlsx")

rename <- dplyr::rename

#renaming bulls eye items from Swedish to English
BullsEye_merged_data_final_without_comp_level <-
  BullsEye_merged_data_final_without_comp_level %>%
  rename(
    BE_Competition = BE_Tavling,
    BE_Training = BE_Traning,
    BE_PreparationRecovery = BE_Forberedelser,
    BE_Life = BE_Livet,
    BE_Obstacles_Sport = BE_Hinder_Idrotten,
    BE_Obstacles_Life = BE_Hinder_Livet,
    BE_Competition_test2 = BE_Tavling_test2,
    BE_Training_test2 = BE_Traning_test2,
    BE_PreparationRecovery_test2 = BE_Forberedelser_test2,
    BE_Life_test2 = BE_Livet_test2,
    BE_Obstacles_Sport_test2 = BE_Hinder_Idrotten_test2,
    BE_Obstacles_Life_test2 = BE_Hinder_Livet_test2
  )

#adding competition level
comp_level_BE <- read_excel("comp_level_bulls eye.xlsx")

BullsEye_merged_data_final <- full_join(BullsEye_merged_data_final_without_comp_level,
                                        comp_level_BE,
                                        by = "Användarnamn")

#data corrections for BullsEye_merged_data_final
#age reported as year of birth
BullsEye_merged_data_final$Demografi.Alder[BullsEye_merged_data_final$Användarnamn == "1054fmsc" &
                                             BullsEye_merged_data_final$Demografi.Alder == 95] <- 24
BullsEye_merged_data_final$Demografi.Alder[BullsEye_merged_data_final$Användarnamn == "1180smpq" &
                                             BullsEye_merged_data_final$Demografi.Alder == 91] <- 28
BullsEye_merged_data_final$Demografi.Alder[BullsEye_merged_data_final$Användarnamn == "1183qhst" &
                                             BullsEye_merged_data_final$Demografi.Alder == 93] <- 26

#participant stopped practice sport by test2 so values at test2 therefore not included
variables_to_na_1086fcfd <- c(
  "BE_Competition_test2",
  "BE_Training_test2",
  "BE_PreparationRecovery_test2",
  "BE_Life_test2",
  "BE_Obstacles_Sport_test2",
  "BE_Obstacles_Life_test2",
  "SMS-2.1_test2",
  "SMS-2.2_test2",
  "SMS-2.3_test2",
  "SMS-2.4_test2",
  "SMS-2.5_test2",
  "SMS-2.6_test2",
  "SMS-2.7_test2",
  "SMS-2.8_test2",
  "SMS-2.9_test2",
  "SMS-2.10_test2",
  "SMS-2.11_test2",
  "SMS-2.12_test2",
  "SMS-2.13_test2",
  "SMS-2.14_test2",
  "SMS-2.15_test2",
  "SMS-2.16_test2",
  "SMS-2.17_test2",
  "SMS-2.18_test2"
)

BullsEye_merged_data_final[BullsEye_merged_data_final$Användarnamn == "1086fcfd", variables_to_na_1086fcfd] <- NA

#rename Valuing Questionnaire items abbrevation from VLQ to VQ
BullsEye_merged_data_final <- BullsEye_merged_data_final %>%
  rename_with( ~ gsub("^VLQ", "VQ", .), starts_with("VLQ"))


#create date and time difference variables, for survey 1 and survey 2
# Ensure the variables are in POSIXct format
BullsEye_merged_data_final$Färdig_surveys_test1 <- as.POSIXct(BullsEye_merged_data_final$Färdig_surveys_test1, format = "%Y-%m-%d %H:%M")
BullsEye_merged_data_final$Färdig_surveys_test2 <- as.POSIXct(BullsEye_merged_data_final$Färdig_surveys_test2, format = "%Y-%m-%d %H:%M")

# Calculate the time difference in hours (for survey 1 and survey 2)
BullsEye_merged_data_final$time_latency_survey_1_2_hours <- as.numeric(
  difftime(
    BullsEye_merged_data_final$Färdig_surveys_test2,
    BullsEye_merged_data_final$Färdig_surveys_test1,
    units = "hours"
  )
)

# Calculate the time difference in days (for survey 1 and survey 2)
BullsEye_merged_data_final$time_latency_survey_1_2_days <- as.numeric(
  difftime(
    BullsEye_merged_data_final$Färdig_surveys_test2,
    BullsEye_merged_data_final$Färdig_surveys_test1,
    units = "days"
  )
)



#create date and time difference variables, for survey 1 and bulls eye 1
# Ensure the variables are in POSIXct format
BullsEye_merged_data_final$Färdig_surveys_test1 <- as.POSIXct(BullsEye_merged_data_final$Färdig_surveys_test1, format = "%Y-%m-%d %H:%M")
BullsEye_merged_data_final$Skapat <- as.POSIXct(BullsEye_merged_data_final$Skapat, format = "%Y-%m-%d %H:%M:%S")

# Calculate the time difference in hours (for survey 1 and bulls eye 1)
BullsEye_merged_data_final$time_latency_survey1_bulls_eye1_hours <- as.numeric(
  difftime(
    BullsEye_merged_data_final$Skapat,
    BullsEye_merged_data_final$Färdig_surveys_test1,
    units = "hours"
  )
)

# Calculate the time difference in days (for survey 1 and bulls eye 1)
BullsEye_merged_data_final$time_latency_survey1_bulls_eye1_days <- as.numeric(
  difftime(
    BullsEye_merged_data_final$Skapat,
    BullsEye_merged_data_final$Färdig_surveys_test1,
    units = "days"
  )
)


#create date and time difference variables, for bulls eye 1 and bulls eye 2
# Ensure the variables are in POSIXct format
BullsEye_merged_data_final$Skapat <- as.POSIXct(BullsEye_merged_data_final$Skapat, format = "%Y-%m-%d %H:%M:%S")
BullsEye_merged_data_final$Skapat_test2 <- as.POSIXct(BullsEye_merged_data_final$Skapat_test2, format = "%Y-%m-%d %H:%M:%S")

# Calculate the time difference in hours (for bulls eye 1 and bulls eye 2)
BullsEye_merged_data_final$time_latency_bulls_eye1_2_hours <- as.numeric(
  difftime(
    BullsEye_merged_data_final$Skapat_test2,
    BullsEye_merged_data_final$Skapat,
    units = "hours"
  )
)

# Calculate the time difference in days (for bulls eye 1 and bulls eye 2)
BullsEye_merged_data_final$time_latency_bulls_eye1_2_days <- as.numeric(
  difftime(
    BullsEye_merged_data_final$Skapat_test2,
    BullsEye_merged_data_final$Skapat,
    units = "days"
  )
)


#have checked and all date and time variables are correct 

filter <- dplyr::filter

#hockeystudien bulls eye data:
sheet1 <- read_excel("hockeystudien intervention_bulls eye.xlsx", sheet = 1)
sheet1_filtered <- sheet1 %>%
  filter(`Instansens namn` == "Vecka 1") %>%
  distinct(`Användarnamn`, .keep_all = TRUE)

sheet2 <- read_excel("hockeystudien intervention_bulls eye.xlsx", sheet = 2)
sheet2_filtered <- sheet2 %>%
  filter(`Instansens namn` == "Vecka 1") %>%
  distinct(`Användarnamn`, .keep_all = TRUE)

merged_data_hockeystudien_BE <- full_join(sheet1_filtered, sheet2_filtered, by = "Användarnamn")

#hockeystudien demographics
hockeystudien_intervention_demographics <- read_excel("hockeystudien intervention_demographics.xlsx", sheet = 1)

#corecting age for a participant in the df
hockeystudien_intervention_demographics$Demografi.Alder[hockeystudien_intervention_demographics$Användarnamn == "1077dhrv" &
                                                          hockeystudien_intervention_demographics$Demografi.Alder == 6] <- 26


#merge df BE and demographics from hockeystudien
merged_data_hockeystudien_BE_demographics <- full_join(merged_data_hockeystudien_BE,
                                                       hockeystudien_intervention_demographics,
                                                       by = "Användarnamn")

#following in the hockeystudien study already exists in the bulls eye study, therefore mask:
#mask 1029spns, 1061bgmb, 1043gsnq, 1046pkbz, 1062xhqd, 1070mfqq

usernames_to_exclude <- c("1029spns",
                          "1061bgmb",
                          "1043gsnq",
                          "1046pkbz",
                          "1062xhqd",
                          "1070mfqq")

filtered_data <- merged_data_hockeystudien_BE_demographics %>%
  filter(!(`Användarnamn` %in% usernames_to_exclude))


hockeystudien_bulls_eye_modified_data <- filtered_data %>%
  mutate(`Användarnamn` = paste0(`Användarnamn`, "_hockeystudien")) %>%  # Append "_hockeystudien" to each value
  rename(
    BE_Competition = Tavling,
    BE_Training = Traning,
    BE_PreparationRecovery = Forberedelser,
    BE_Life = Livet,
    BE_Obstacles_Sport = HinderSkala1,
    BE_Obstacles_Life = HinderSkala2
  )




#preparations for merging datasets into superfinal version (bulls eye + hockeystudien_bulls eye)

hockeystudien_bulls_eye_modified_data <- hockeystudien_bulls_eye_modified_data %>%
  mutate(from_hockeystudien = 1)

BullsEye_merged_data_final <- BullsEye_merged_data_final %>%
  mutate(from_hockeystudien = 0)

bulls_eye_merged_superfinal_data <- bind_rows(hockeystudien_bulls_eye_modified_data,
                                              BullsEye_merged_data_final) #use for analyses

#BE items (dartboard), adjust (-1) to enable Rasch analysis
BE_adjust <- c(
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Competition_test2",
  "BE_Training_test2",
  "BE_PreparationRecovery_test2",
  "BE_Life_test2"
)

bulls_eye_merged_superfinal_data[, BE_adjust] <- lapply(bulls_eye_merged_superfinal_data[, BE_adjust], function(x)
  x - 1)

#"Hinder", reverse items (low value = low in values-based living and compatible with other values items) and adjust (-1) to enable Rasch analysis
Hinder_reverse_and_adjust <- c(
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life",
  "BE_Obstacles_Sport_test2",
  "BE_Obstacles_Life_test2"
)

bulls_eye_merged_superfinal_data[, Hinder_reverse_and_adjust] <- lapply(bulls_eye_merged_superfinal_data[, Hinder_reverse_and_adjust], function(x)
  8 - x - 1)


#ADAQ, reverse (higher value = higher PF) and adjust (-1) to enable Rasch analysis
ADAQ_reverse_and_adjust <- c("ADAQ.1",
                             "ADAQ.2",
                             "ADAQ.3",
                             "ADAQ.4",
                             "ADAQ.5",
                             "ADAQ.6",
                             "ADAQ.7")

bulls_eye_merged_superfinal_data[, ADAQ_reverse_and_adjust] <- lapply(bulls_eye_merged_superfinal_data[, ADAQ_reverse_and_adjust], function(x)
  8 - x - 1)


#SWLS, adjust (-1) to enable Rasch analysis
SWLS_adjust <- c("SWLS.1", "SWLS.2", "SWLS.3", "SWLS.4", "SWLS.5")

bulls_eye_merged_superfinal_data[, SWLS_adjust] <- lapply(bulls_eye_merged_superfinal_data[, SWLS_adjust], function(x)
  x - 1)

#VQ scale already starts at zero and does therefore not need to be adjusted for Rasch analysis

2 Bulls-Eye for Athletes: Rasch analysis

Code
Code
#creating a data frame with the bulls eye items
BE_6items_data_rasch <- bulls_eye_merged_superfinal_data[, c(
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life"
)]
Code
#investigate missing data
RImissing(BE_6items_data_rasch)

Code
#the following means that participants have answered at least one item
filter <- dplyr::filter

min.responses <- 1

BE_6items_data_rasch <- BE_6items_data_rasch %>%
  filter(6 - rowSums(is.na(
    select(
      .,
      "BE_Competition",
      "BE_Training",
      "BE_PreparationRecovery",
      "BE_Life",
      "BE_Obstacles_Sport",
      "BE_Obstacles_Life"
    )
  )) >= min.responses)
Code
#tile plot (response category frequency for the items)
RItileplot(BE_6items_data_rasch)

Code
#investigating reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI
RItif(na.omit(BE_6items_data_rasch))

Code
RItif(na.omit(BE_6items_data_rasch), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(BE_6items_data_rasch, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue sig
overall 44.7 35 0.13
Code
#conditional item fit (assessing unidimensionaliity)
simfit1_BE_6items <- RIgetfit(BE_6items_data_rasch, iterations = 100, cpu = 10) 

RIitemfit(BE_6items_data_rasch, simfit1_BE_6items, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
BE_Competition 1.146 [0.823, 1.254] 1.181 [0.781, 1.366] no misfit no misfit -0.18
BE_Training 1.009 [0.817, 1.201] 1.007 [0.805, 1.231] no misfit no misfit -0.20
BE_PreparationRecovery 0.874 [0.757, 1.246] 0.883 [0.773, 1.262] no misfit no misfit -0.19
BE_Life 1.014 [0.831, 1.212] 1.016 [0.83, 1.271] no misfit no misfit -0.37
BE_Obstacles_Sport 0.972 [0.78, 1.161] 0.948 [0.792, 1.164] no misfit no misfit 0.13
BE_Obstacles_Life 1 [0.782, 1.181] 0.994 [0.776, 1.172] no misfit no misfit -0.42
Note:
MSQ values based on conditional calculations (n = 146 complete cases).
Simulation based thresholds from 85 simulated datasets.
Code
#conditional item carachteristic curves 
ICCplot(
  as.data.frame(BE_6items_data_rasch),
  itemnumber = 1:4,
  itemdescrip = c("Competition", "Training", "Preparation & Recovery", "Life"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(BE_6items_data_rasch),
  itemnumber = 5:6,
  itemdescrip = c("Obstacles Sport", "Obstacles Life"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
#PCA of residuals (assessing unidimensionality)
RIpcmPCA(na.omit(BE_6items_data_rasch))
Eigenvalues Proportion of variance
1.86 32.1%
1.38 22.5%
1.07 16.5%
0.93 16.1%
0.75 12.6%
Code
#item restscore
RIrestscore(BE_6items_data_rasch)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
BE_Competition 0.33 0.33 0.00 0.992 0.02 -0.18
BE_Training 0.29 0.29 0.00 0.992 0.01 -0.20
BE_PreparationRecovery 0.38 0.30 0.08 0.992 0.02 -0.19
BE_Life 0.30 0.31 0.01 0.992 -0.16 -0.37
BE_Obstacles_Sport 0.33 0.31 0.02 0.992 0.33 0.13
BE_Obstacles_Life 0.32 0.32 0.00 0.992 -0.22 -0.42
Code
#loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(na.omit(BE_6items_data_rasch))

Code
#residual correlations
resid_cor_BE_6items <- RIgetResidCor(BE_6items_data_rasch, iterations = 800, cpu = 10)

RIresidcorr(BE_6items_data_rasch, cutoff = resid_cor_BE_6items$p99)
BE_Competition BE_Training BE_PreparationRecovery BE_Life BE_Obstacles_Sport BE_Obstacles_Life
BE_Competition
BE_Training 0.09
BE_PreparationRecovery -0.08 -0.11
BE_Life -0.22 -0.11 -0.16
BE_Obstacles_Sport -0.23 -0.25 0 -0.22
BE_Obstacles_Life -0.44 -0.26 -0.09 0.04 0.12
Note:
Relative cut-off value is 0.145, which is 0.274 above the average correlation (-0.129).
Correlations above the cut-off are highlighted in red text.
Code
#partial gamma coefficients (assessing local independence)

clean_names <- janitor::clean_names
filter <- dplyr::filter

RIpartgamLD(BE_6items_data_rasch) #shows only significant results
Item 1 Item 2 Partial gamma SE Lower CI Upper CI Adjusted p-value (BH)
Code
partgam_LD(as.data.frame(BE_6items_data_rasch)) #shows all output, also non-significant results
                    Item1                  Item2   gamma     se pvalue padj.BH
1          BE_Competition            BE_Training  0.2089 0.1315 0.1121  1.0000
2          BE_Competition BE_PreparationRecovery  0.1289 0.1175 0.2727  1.0000
3          BE_Competition                BE_Life  0.0097 0.1262 0.9387  1.0000
4          BE_Competition     BE_Obstacles_Sport -0.0534 0.1332 0.6882  1.0000
5          BE_Competition      BE_Obstacles_Life -0.2618 0.1222 0.0322  0.9663
6             BE_Training BE_PreparationRecovery -0.0806 0.1123 0.4727  1.0000
7             BE_Training                BE_Life -0.0178 0.1293 0.8905  1.0000
8             BE_Training     BE_Obstacles_Sport -0.1873 0.1265 0.1387  1.0000
9             BE_Training      BE_Obstacles_Life -0.1478 0.1176 0.2089  1.0000
10 BE_PreparationRecovery                BE_Life -0.0254 0.1281 0.8429  1.0000
11 BE_PreparationRecovery     BE_Obstacles_Sport  0.0025 0.1223 0.9839  1.0000
12 BE_PreparationRecovery      BE_Obstacles_Life -0.0724 0.1225 0.5546  1.0000
13                BE_Life     BE_Obstacles_Sport -0.1949 0.1362 0.1525  1.0000
14                BE_Life      BE_Obstacles_Life  0.1802 0.1195 0.1316  1.0000
15     BE_Obstacles_Sport      BE_Obstacles_Life  0.2875 0.1141 0.0117  0.3524
   sig   lower   upper
1      -0.0488  0.4666
2      -0.1014  0.3592
3      -0.2377  0.2571
4      -0.3144  0.2075
5      -0.5013 -0.0222
6      -0.3008  0.1395
7      -0.2713  0.2357
8      -0.4353  0.0606
9      -0.3784  0.0828
10     -0.2764  0.2257
11     -0.2373  0.2422
12     -0.3125  0.1677
13     -0.4620  0.0721
14     -0.0540  0.4145
15      0.0639  0.5112

                    Item1                  Item2   gamma     se pvalue padj.BH
1             BE_Training         BE_Competition  0.1952 0.1150 0.0897  1.0000
2  BE_PreparationRecovery         BE_Competition  0.1799 0.1093 0.0998  1.0000
3  BE_PreparationRecovery            BE_Training  0.0317 0.1119 0.7769  1.0000
4                 BE_Life         BE_Competition  0.0769 0.1051 0.4642  1.0000
5                 BE_Life            BE_Training  0.0201 0.1171 0.8639  1.0000
6                 BE_Life BE_PreparationRecovery -0.1082 0.1190 0.3629  1.0000
7      BE_Obstacles_Sport         BE_Competition -0.1304 0.1179 0.2684  1.0000
8      BE_Obstacles_Sport            BE_Training -0.1291 0.1367 0.3450  1.0000
9      BE_Obstacles_Sport BE_PreparationRecovery  0.1530 0.1321 0.2467  1.0000
10     BE_Obstacles_Sport                BE_Life -0.1940 0.1234 0.1161  1.0000
11      BE_Obstacles_Life         BE_Competition -0.3671 0.1071 0.0006  0.0182
12      BE_Obstacles_Life            BE_Training -0.1755 0.1137 0.1228  1.0000
13      BE_Obstacles_Life BE_PreparationRecovery -0.0439 0.1248 0.7248  1.0000
14      BE_Obstacles_Life                BE_Life  0.1898 0.1130 0.0930  1.0000
15      BE_Obstacles_Life     BE_Obstacles_Sport  0.3054 0.1047 0.0035  0.1061
   sig   lower   upper
1      -0.0302  0.4206
2      -0.0343  0.3941
3      -0.1876  0.2510
4      -0.1291  0.2829
5      -0.2094  0.2495
6      -0.3414  0.1250
7      -0.3614  0.1006
8      -0.3971  0.1389
9      -0.1059  0.4119
10     -0.4359  0.0480
11   * -0.5769 -0.1573
12     -0.3984  0.0474
13     -0.2885  0.2006
14     -0.0316  0.4113
15      0.1002  0.5107
Code
#targeting
RItargeting(BE_6items_data_rasch)

Code
#item hierarchy
itemlabels <- data.frame(itemnr = names(BE_6items_data_rasch), item = NA)

RIitemHierarchy(BE_6items_data_rasch)

Code
#threshold locations 
RIitemparams(BE_6items_data_rasch)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Threshold 6 Item location
BE_Competition 0.15 -0.58 -0.92 0.45 0.72 0.30 0.02
BE_Training -0.58 -1.91 -0.71 0.55 0.54 2.13 0.01
BE_PreparationRecovery -1.56 -0.84 -0.50 0.21 0.62 2.17 0.02
BE_Life -1.37 -1.24 -0.17 0.09 0.83 0.90 -0.16
BE_Obstacles_Sport -1.02 -1.16 0.40 0.38 1.11 2.29 0.33
BE_Obstacles_Life -2.03 -1.18 0.16 0.42 0.61 0.73 -0.22
Note:
Item location is the average of the thresholds for each item.
Code
#analysis of response categories
RIitemCats(BE_6items_data_rasch, items = "all", xlims = c(-5, 5), legend = "left")

Based on the analysis ofResponse categories will

2.1 Recoding response categories for Bulls-Eye for Athletes

Based on the analysis of the response categories and threshold locations, scale steps were merged and recoded so that each kept scale step had the highest response probability at some point on the person location scale, and also so that all threshold locations were in the right order for all items. Recoded items were given new names with the abbrevation “_recoded”.

Code
#adding recoded bulls eye items
recode <- dplyr::recode

bulls_eye_merged_superfinal_data <- bulls_eye_merged_superfinal_data %>%
  mutate(
    # Recoding BE_Competition (0+1; 2+3; 4+5)
    BE_Competition_recoded = recode(
      BE_Competition,
      `0` = 0, `1` = 0,
      `2` = 1, `3` = 1,
      `4` = 2, `5` = 2,
      `6` = 3
    ),
    # Recoding BE_Training (0+1; 3+4)
    BE_Training_recoded = recode(
      BE_Training,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2, `4` = 2,
      `5` = 3,
      `6` = 4
    ),
    #Recoding BE_PreparationRecovery(0+1)
    BE_PreparationRecovery_recoded = recode(
      BE_PreparationRecovery,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2,
      `4` = 3,
      `5` = 4,
      `6` = 5
    ),
    #Recoding BE_Life (0+1)
    BE_Life_recoded = recode(
      BE_Life,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2,
      `4` = 3,
      `5` = 4,
      `6` = 5
    ),
    #Recoding BE_Obstacles_Sport (0+1; 3+4; 5+6)
    BE_Obstacles_Sport_recoded = recode(
      BE_Obstacles_Sport,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2, `4` = 2,
      `5` = 3, `6` = 3
    ),
    #Recoding BE_Obstacles_Life (0+1)
    BE_Obstacles_Life_recoded = recode(
      BE_Obstacles_Life,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2,
      `4` = 3,
      `5` = 4,
      `6` = 5
    )
  )

2.2 Re-analysis of Bulls Eye recoded version

Code
#preparing item data frame for BE recoded

BE_6items_data_rasch_recoded <- bulls_eye_merged_superfinal_data[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

min.responses <- 1

BE_6items_data_rasch_recoded <- BE_6items_data_rasch_recoded %>%
  filter(6 - rowSums(is.na(
    select(
      .,
      "BE_Competition_recoded",
      "BE_Training_recoded",
      "BE_PreparationRecovery_recoded",
      "BE_Life_recoded",
      "BE_Obstacles_Sport_recoded",
      "BE_Obstacles_Life_recoded"
    )
  )) >= min.responses)
Code
#re-analysis of response categories using BE recoded
RIitemCats(BE_6items_data_rasch_recoded, items = "all", xlims = c(-5, 5), legend = "left")

Code
#tile plot (response category frequency for the recoded items)
RItileplot(BE_6items_data_rasch_recoded)

Code
#re-analysis reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI
RItif(na.omit(BE_6items_data_rasch_recoded))

Code
RItif(na.omit(BE_6items_data_rasch_recoded), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(BE_6items_data_rasch_recoded, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue  sig 
overall 57.9 24 0.00013  ***
Code
#re-analysis of conditional item fit (assessing unidimensionaliity)
simfit1_BE_6items_recoded <- RIgetfit(BE_6items_data_rasch_recoded, iterations = 100, cpu = 10) 

RIitemfit(BE_6items_data_rasch_recoded, simfit1_BE_6items_recoded, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
BE_Competition_recoded 0.97 [0.804, 1.125] 0.966 [0.795, 1.183] no misfit no misfit -0.05
BE_Training_recoded 1.042 [0.784, 1.22] 1.053 [0.784, 1.273] no misfit no misfit 0.03
BE_PreparationRecovery_recoded 0.959 [0.774, 1.211] 0.965 [0.798, 1.238] no misfit no misfit 0.18
BE_Life_recoded 1.076 [0.811, 1.156] 1.078 [0.81, 1.244] no misfit no misfit -0.07
BE_Obstacles_Sport_recoded 0.936 [0.801, 1.205] 0.967 [0.792, 1.241] no misfit no misfit 0.00
BE_Obstacles_Life_recoded 1.005 [0.795, 1.221] 1.021 [0.778, 1.223] no misfit no misfit -0.02
Note:
MSQ values based on conditional calculations (n = 146 complete cases).
Simulation based thresholds from 100 simulated datasets.
Code
#re-analysis of conditional item carachteristic curves
ICCplot(
  as.data.frame(BE_6items_data_rasch_recoded),
  itemnumber = 1:4,
  itemdescrip = c(
    "Competition_recoded",
    "Training_recoded",
    "Preparation & Recovery_recoded",
    "Life_recoded"
  ),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(BE_6items_data_rasch_recoded),
  itemnumber = 5:6,
  itemdescrip = c("Obstacles Sport_recoded", "Obstacles Life_recoded"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
#re-analyzing PCA of residuals (assessing unidimensionality)
RIpcmPCA(na.omit(BE_6items_data_rasch_recoded))
Eigenvalues Proportion of variance
1.80 30.2%
1.45 24.1%
1.08 17.4%
0.97 16.5%
0.68 11.6%
Code
#re-analyzing item restscores
RIrestscore(BE_6items_data_rasch_recoded)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
BE_Competition_recoded 0.37 0.31 0.06 0.86 -0.07 -0.05
BE_Training_recoded 0.29 0.31 0.02 0.86 0.01 0.03
BE_PreparationRecovery_recoded 0.36 0.35 0.01 0.86 0.16 0.18
BE_Life_recoded 0.31 0.36 0.05 0.86 -0.09 -0.07
BE_Obstacles_Sport_recoded 0.40 0.30 0.10 0.86 -0.02 0.00
BE_Obstacles_Life_recoded 0.34 0.37 0.03 0.86 -0.04 -0.02
Code
#re-analyzing loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(na.omit(BE_6items_data_rasch_recoded))

Code
#re-analysis of residual correlations
resid_cor_BE_6items_recoded <- RIgetResidCor(BE_6items_data_rasch_recoded,
                                             iterations = 800,
                                             cpu = 10)

RIresidcorr(BE_6items_data_rasch_recoded, cutoff = resid_cor_BE_6items_recoded$p99)
BE_Competition_recoded BE_Training_recoded BE_PreparationRecovery_recoded BE_Life_recoded BE_Obstacles_Sport_recoded BE_Obstacles_Life_recoded
BE_Competition_recoded
BE_Training_recoded 0.2
BE_PreparationRecovery_recoded 0 -0.08
BE_Life_recoded -0.14 -0.14 -0.27
BE_Obstacles_Sport_recoded -0.08 -0.26 -0.04 -0.23
BE_Obstacles_Life_recoded -0.37 -0.25 -0.23 -0.11 0.11
Note:
Relative cut-off value is 0.152, which is 0.278 above the average correlation (-0.126).
Correlations above the cut-off are highlighted in red text.
Code
#re-analysis of partial gamma coefficients (assessing local independence)
RIpartgamLD(BE_6items_data_rasch_recoded) #shows only significant results
Item 1 Item 2 Partial gamma SE Lower CI Upper CI Adjusted p-value (BH)
BE_Competition_recoded BE_Training_recoded 0.444 0.141 0.169 0.720 0.047
BE_Training_recoded BE_Competition_recoded 0.419 0.123 0.177 0.661 0.020
BE_Obstacles_Sport_recoded BE_Obstacles_Life_recoded 0.391 0.111 0.174 0.609 0.012
Code
partgam_LD(as.data.frame(BE_6items_data_rasch_recoded)) #shows all output, also non-significant results
                            Item1                          Item2   gamma     se
1          BE_Competition_recoded            BE_Training_recoded  0.4444 0.1407
2          BE_Competition_recoded BE_PreparationRecovery_recoded  0.1895 0.1266
3          BE_Competition_recoded                BE_Life_recoded  0.0435 0.1249
4          BE_Competition_recoded     BE_Obstacles_Sport_recoded -0.0163 0.1471
5          BE_Competition_recoded      BE_Obstacles_Life_recoded -0.2319 0.1338
6             BE_Training_recoded BE_PreparationRecovery_recoded  0.0354 0.1316
7             BE_Training_recoded                BE_Life_recoded  0.0000 0.1391
8             BE_Training_recoded     BE_Obstacles_Sport_recoded -0.1847 0.1606
9             BE_Training_recoded      BE_Obstacles_Life_recoded -0.1848 0.1267
10 BE_PreparationRecovery_recoded                BE_Life_recoded -0.0282 0.1235
11 BE_PreparationRecovery_recoded     BE_Obstacles_Sport_recoded -0.0559 0.1475
12 BE_PreparationRecovery_recoded      BE_Obstacles_Life_recoded -0.1255 0.1218
13                BE_Life_recoded     BE_Obstacles_Sport_recoded -0.2246 0.1427
14                BE_Life_recoded      BE_Obstacles_Life_recoded  0.0979 0.1189
15     BE_Obstacles_Sport_recoded      BE_Obstacles_Life_recoded  0.3915 0.1107
   pvalue padj.BH sig   lower  upper
1  0.0016  0.0474   *  0.1688 0.7201
2  0.1344  1.0000     -0.0586 0.4377
3  0.7278  1.0000     -0.2014 0.2883
4  0.9119  1.0000     -0.3047 0.2721
5  0.0831  1.0000     -0.4941 0.0304
6  0.7882  1.0000     -0.2226 0.2933
7  1.0000  1.0000     -0.2727 0.2727
8  0.2503  1.0000     -0.4995 0.1301
9  0.1449  1.0000     -0.4332 0.0636
10 0.8193  1.0000     -0.2702 0.2138
11 0.7049  1.0000     -0.3450 0.2333
12 0.3027  1.0000     -0.3643 0.1132
13 0.1155  1.0000     -0.5042 0.0551
14 0.4103  1.0000     -0.1351 0.3309
15 0.0004  0.0122   *  0.1745 0.6085

                            Item1                          Item2   gamma     se
1             BE_Training_recoded         BE_Competition_recoded  0.4190 0.1233
2  BE_PreparationRecovery_recoded         BE_Competition_recoded  0.1256 0.1294
3  BE_PreparationRecovery_recoded            BE_Training_recoded  0.0671 0.1341
4                 BE_Life_recoded         BE_Competition_recoded  0.0246 0.1194
5                 BE_Life_recoded            BE_Training_recoded -0.0122 0.1206
6                 BE_Life_recoded BE_PreparationRecovery_recoded -0.1278 0.1234
7      BE_Obstacles_Sport_recoded         BE_Competition_recoded -0.0646 0.1472
8      BE_Obstacles_Sport_recoded            BE_Training_recoded -0.1852 0.1605
9      BE_Obstacles_Sport_recoded BE_PreparationRecovery_recoded  0.0823 0.1416
10     BE_Obstacles_Sport_recoded                BE_Life_recoded -0.1557 0.1396
11      BE_Obstacles_Life_recoded         BE_Competition_recoded -0.4621 0.1084
12      BE_Obstacles_Life_recoded            BE_Training_recoded -0.2538 0.1205
13      BE_Obstacles_Life_recoded BE_PreparationRecovery_recoded -0.1522 0.1196
14      BE_Obstacles_Life_recoded                BE_Life_recoded  0.0776 0.1061
15      BE_Obstacles_Life_recoded     BE_Obstacles_Sport_recoded  0.2752 0.1233
   pvalue padj.BH  sig   lower   upper
1  0.0007  0.0204    *  0.1773  0.6607
2  0.3317  1.0000      -0.1280  0.3793
3  0.6168  1.0000      -0.1957  0.3298
4  0.8365  1.0000      -0.2093  0.2586
5  0.9195  1.0000      -0.2486  0.2242
6  0.3003  1.0000      -0.3697  0.1140
7  0.6608  1.0000      -0.3532  0.2240
8  0.2485  1.0000      -0.4997  0.1293
9  0.5611  1.0000      -0.1952  0.3598
10 0.2648  1.0000      -0.4293  0.1180
11 0.0000  0.0006  *** -0.6747 -0.2496
12 0.0351  1.0000      -0.4900 -0.0177
13 0.2032  1.0000      -0.3867  0.0822
14 0.4648  1.0000      -0.1304  0.2855
15 0.0256  0.7677       0.0336  0.5168
Code
#targeting
RItargeting(BE_6items_data_rasch_recoded)

Code
#item hierarchy
itemlabels <- data.frame(itemnr = names(BE_6items_data_rasch_recoded), item = NA)

RIitemHierarchy(BE_6items_data_rasch_recoded)

Code
#threshold locations 
RIitemparams(BE_6items_data_rasch_recoded)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Item location
BE_Competition_recoded -1.20 0.04 0.96 NA NA -0.07
BE_Training_recoded -1.66 -1.44 1.17 1.98 NA 0.01
BE_PreparationRecovery_recoded -0.96 -0.78 -0.01 0.47 2.08 0.16
BE_Life_recoded -1.31 -0.46 -0.13 0.67 0.81 -0.09
BE_Obstacles_Sport_recoded -1.09 -0.43 1.45 NA NA -0.02
BE_Obstacles_Life_recoded -1.37 -0.11 0.21 0.46 0.63 -0.04
Note:
Item location is the average of the thresholds for each item.

2.3 Differential Item Functioning (DIF) and test-retest investigation

DIF investigates invariance, and more specifically whether item functioning varies between subgroups of participants or not. DIF variables to be investigated are age, sex, competition level, sport type (individual/team sport), and high/low scorers on the instrument. All DIF variables will be transformed into binary categories. Test-retest reliability will also be investigated within the framework of DIF. Below is data preparation and analysis for each DIF variable.

2.3.1 DIF analyses

Code
#Age DIF variable data preparation

median_age <- median(bulls_eye_merged_superfinal_data$Demografi.Alder)
print(median_age)
[1] 23
Code
#median age is 23 

specific_value <- 23

count_value_median_age <- sum(bulls_eye_merged_superfinal_data$Demografi.Alder == specific_value, na.rm = TRUE)

print(count_value_median_age)
[1] 14
Code
#14 participants is 23 years of age 

range(bulls_eye_merged_superfinal_data$Demografi.Alder)
[1] 15 50
Code
#age range is 15 to 50 years

#age variable is split by median using an inclusive split in which
#participants aged 23 is included in the younger group
#creating age_group variable

bulls_eye_merged_superfinal_data$age_group <- ifelse(bulls_eye_merged_superfinal_data$Demografi.Alder <= median_age,
                                                     "Younger",
                                                     "Older")

dif.age <- factor(bulls_eye_merged_superfinal_data$age_group)

#Age DIF analysis
RIpartgamDIF(BE_6items_data_rasch_recoded, dif.age) 
[1] "No statistically significant DIF found."
Code
#Sex DIF variable data preparation
table(bulls_eye_merged_superfinal_data$Demografi.Kon)

 1  2 
66 90
Code
#66 females and 90 males 

dif.sex <- factor(bulls_eye_merged_superfinal_data$Demografi.Kon)

#Sex DIF analysis
RIpartgamDIF(BE_6items_data_rasch_recoded, dif.sex)
[1] "No statistically significant DIF found."
Code
#creating binary variable for competition level

#available competition levels in the current sample are
#junior elite, sectional, national, international
#these are categorized as "Lower level competition" (junior elite + sectional)
#and "Higher competition level" (national + international) 

bulls_eye_merged_superfinal_data$comp_level_binary <- ifelse(
  bulls_eye_merged_superfinal_data$comp_level %in% c("junior elite", "sectional"), 
  "Lower competition level", 
  "Higher competition level"
)

table(bulls_eye_merged_superfinal_data$comp_level_binary)

Higher competition level  Lower competition level 
                      93                       63
Code
#93 higher comp level, 63 lower comp level

dif.comp_level <- factor(bulls_eye_merged_superfinal_data$comp_level_binary)

#Competition level DIF analysis
RIpartgamDIF(BE_6items_data_rasch_recoded, dif.comp_level) 
[1] "No statistically significant DIF found."
Code
#Sport type DIF variable data preparation (individual vs team sport)

table(bulls_eye_merged_superfinal_data$Demografi.Idrotttyp)

  1   2 
 19 137
Code
#19 individual and 137 team sport athletes

dif.sport_type <- factor(bulls_eye_merged_superfinal_data$Demografi.Idrotttyp)

#Sport type DIF analysis
RIpartgamDIF(BE_6items_data_rasch_recoded, dif.sport_type) 
[1] "No statistically significant DIF found."
Code
#DIF analysis for high vs low scorers on the bulls eye instrument (group 1 = low scorers; group 2 = high scorers)
iarm::item_obsexp(PCM(BE_6items_data_rasch_recoded))
Score group 1: 
                               mean obs mean exp std.res sig
BE_Competition_recoded          1.127    1.165   -0.427     
BE_Training_recoded             1.810    1.743    0.855     
BE_PreparationRecovery_recoded  1.873    1.803    0.606     
BE_Life_recoded                 1.899    1.792    0.908     
BE_Obstacles_Sport_recoded      1.190    1.261   -0.808     
BE_Obstacles_Life_recoded       1.443    1.302    1.386     

Score group 2: 
                               mean obs mean exp std.res sig
BE_Competition_recoded          2.000    1.955    0.463     
BE_Training_recoded             2.239    2.318   -0.930     
BE_PreparationRecovery_recoded  3.164    3.216   -0.446     
BE_Life_recoded                 3.239    3.364   -0.990     
BE_Obstacles_Sport_recoded      2.030    1.946    0.954     
BE_Obstacles_Life_recoded       3.448    3.320    0.937

2.3.2 Test-retest reliability of Bulls-Eye for Athletes using DIF analysis

Code
#Test-retest reliability using DIF, data preparation and analysis 

#bulls eye "test2" equals time point 2 testing in the code below

rename <- dplyr::rename

#creating data frame with bulls eye test 2 items selected
BE_6items_data_rasch_test2 <- bulls_eye_merged_superfinal_data[, c(
  "Användarnamn",
  "BE_Competition_test2",
  "BE_Training_test2",
  "BE_PreparationRecovery_test2",
  "BE_Life_test2",
  "BE_Obstacles_Sport_test2",
  "BE_Obstacles_Life_test2"
)]

columns_to_check <- c(
  "BE_Competition_test2",
  "BE_Training_test2",
  "BE_PreparationRecovery_test2",
  "BE_Life_test2",
  "BE_Obstacles_Sport_test2",
  "BE_Obstacles_Life_test2"
)

#filter rows with at least one non-missing value in the specified columns
BE_6items_data_rasch_test2_filtered <- BE_6items_data_rasch_test2[rowSums(!is.na(BE_6items_data_rasch_test2[columns_to_check])) > 0, ]
#69 participants with a bulls eye time point 2 test

#rename time point 2 items as time point 1 + create new IDs for time point 2 tests 
BE_6items_data_rasch_test2_filtered_renamed_items <- BE_6items_data_rasch_test2_filtered %>%
  mutate(`Användarnamn` = paste0(`Användarnamn`, "_duplicate")) %>%
  rename(
    "BE_Competition" = "BE_Competition_test2",
    "BE_Training" = "BE_Training_test2",
    "BE_PreparationRecovery" = "BE_PreparationRecovery_test2",
    "BE_Life" = "BE_Life_test2",
    "BE_Obstacles_Sport" = "BE_Obstacles_Sport_test2",
    "BE_Obstacles_Life" = "BE_Obstacles_Life_test2"
  )

#adding DIF variable with value 1 for bulls eye time point 2 (time point 2 tests)
BE_6items_data_rasch_test2_filtered_renamed_items$test2_BE <- 1

#now onto preparing data frame for bulls eye test 1
BE_6items_data_rasch_test1 <- bulls_eye_merged_superfinal_data[, c(
  "Användarnamn",
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life"
)]

#adding DIF variable with value 0 for bulls eye time point 2 (time point 1 tests)
BE_6items_data_rasch_test1$test2_BE <- 0

#merging the two data frames
#in this way time point 2 bulls eye tests can be compared with time point 1
#in which time point 2 tests have been added as "new" individuals in the data (duplicates)
BE_6items_test1_test2_merged <- full_join(BE_6items_data_rasch_test1,
                                          BE_6items_data_rasch_test2_filtered_renamed_items)


#applying the same recoded item structures that came out of the initial BE Rasch analysis (bulls eye recoded version)
BE_6items_test1_test2_merged <- BE_6items_test1_test2_merged %>%
  mutate(
    # Recoding BE_Competition (0+1; 2+3; 4+5)
    BE_Competition_recoded = recode(
      BE_Competition,
      `0` = 0, `1` = 0,
      `2` = 1, `3` = 1,
      `4` = 2, `5` = 2,
      `6` = 3
    ),
    # Recoding BE_Training (0+1; 3+4)
    BE_Training_recoded = recode(
      BE_Training,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2, `4` = 2,
      `5` = 3,
      `6` = 4
    ),
    #Recoding BE_PreparationRecovery(0+1)
    BE_PreparationRecovery_recoded = recode(
      BE_PreparationRecovery,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2,
      `4` = 3,
      `5` = 4,
      `6` = 5
    ),
    #Recoding BE_Life (0+1)
    BE_Life_recoded = recode(
      BE_Life,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2,
      `4` = 3,
      `5` = 4,
      `6` = 5
    ),
    #Recoding BE_Obstacles_Sport (0+1; 3+4; 5+6)
    BE_Obstacles_Sport_recoded = recode(
      BE_Obstacles_Sport,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2, `4` = 2,
      `5` = 3, `6` = 3
    ),
    #Recoding BE_Obstacles_Life (0+1)
    BE_Obstacles_Life_recoded = recode(
      BE_Obstacles_Life,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2,
      `4` = 3,
      `5` = 4,
      `6` = 5
    )
  )


#test-retest DIF investigation
BE_6items_test1_test2_DIF <- BE_6items_test1_test2_merged[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]


dif.test_retest <- factor(BE_6items_test1_test2_merged$test2_BE)

RIpartgamDIF(BE_6items_test1_test2_DIF, dif.test_retest)
[1] "No statistically significant DIF found."

3 Accepance, Defusion and Action Questionnaire-sport (ADAQ): Rasch analysis and BEA validity investigation

The ADAQ instrument is used as part of the validity investigation of BEA. First, ADAQ will undergo Rasch analysis. Then a validity analysis of correlating BEA and ADAQ thetas will be undertaken. Note that ADAQ was not completed by the participants from the intervention study (Hockeystudien), therefore more missing values than in the BEA Rasch analysis are to be expected.

3.1 Rasch analysis of ADAQ

Code
#preparing data frame with ADAQ items 
ADAQ_data_rasch <- bulls_eye_merged_superfinal_data[, c("ADAQ.1",
                                                        "ADAQ.2",
                                                        "ADAQ.3",
                                                        "ADAQ.4",
                                                        "ADAQ.5",
                                                        "ADAQ.6",
                                                        "ADAQ.7")]

#investigate missing values
RImissing(ADAQ_data_rasch)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete ADAQ 
ADAQ_data_rasch_na_omit <- ADAQ_data_rasch %>%
  select(ADAQ.1, ADAQ.2, ADAQ.3, ADAQ.4, ADAQ.5, ADAQ.6, ADAQ.7) %>% 
  na.omit()

RImissing(ADAQ_data_rasch_na_omit)
[1] "No missing data."
Code
#tile plot (response category frequency for the items)
RItileplot(ADAQ_data_rasch_na_omit)

Item ADAQ.2 has no responses in scale step zero, and will therefore be merged with scale step 1.

Code
#merging scale step zero and 1 in ADAQ.2 (creating ADAQ.2_revised) 
bulls_eye_merged_superfinal_data <- bulls_eye_merged_superfinal_data %>%
  mutate(ADAQ.2_revised = recode(
    ADAQ.2,
    `0` = 0, `1` = 0,
    `2` = 1,
    `3` = 2,
    `4` = 3,
    `5` = 4,
    `6` = 5
  ))

#preparing data frame with ADAQ items including ADAQ.2_revised
ADAQ_data_rasch <- bulls_eye_merged_superfinal_data[, c("ADAQ.1",
                                                        "ADAQ.2_revised",
                                                        "ADAQ.3",
                                                        "ADAQ.4",
                                                        "ADAQ.5",
                                                        "ADAQ.6",
                                                        "ADAQ.7")]

#removing missing values 
ADAQ_data_rasch_na_omit <- ADAQ_data_rasch %>%
  select(ADAQ.1, ADAQ.2_revised, ADAQ.3, ADAQ.4, ADAQ.5, ADAQ.6, ADAQ.7) %>% 
  na.omit()

RImissing(ADAQ_data_rasch_na_omit)
[1] "No missing data."
Code
RItileplot(ADAQ_data_rasch_na_omit)

Code
#analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI 
RItif(na.omit(ADAQ_data_rasch_na_omit))

Code
RItif(na.omit(ADAQ_data_rasch_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(ADAQ_data_rasch_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue sig
overall 28.8 40 0.91
Code
#conditional item fit (assessing unidimensionaliity)

simfit1_ADAQ <- RIgetfit(ADAQ_data_rasch_na_omit, iterations = 100, cpu = 10) 

RIitemfit(ADAQ_data_rasch_na_omit, simfit1_ADAQ, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
ADAQ.1 1.254 [0.804, 1.231] 1.296 [0.769, 1.257] 0.023 0.039 -1.22
ADAQ.2_revised 0.737 [0.736, 1.247] 0.751 [0.702, 1.248] no misfit no misfit -1.17
ADAQ.3 0.9 [0.742, 1.239] 0.916 [0.739, 2.351] no misfit no misfit -0.95
ADAQ.4 1.146 [0.735, 1.205] 1.231 [0.731, 1.205] no misfit 0.026 -1.27
ADAQ.5 1.199 [0.775, 1.396] 1.266 [0.77, 1.473] no misfit no misfit -1.14
ADAQ.6 0.745 [0.78, 1.212] 0.73 [0.751, 1.344] 0.035 0.021 -0.87
ADAQ.7 0.975 [0.781, 1.27] 0.972 [0.754, 1.698] no misfit no misfit -1.26
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 59 simulated datasets.
Code
#conditional item charachteristic curves
ICCplot(
  as.data.frame(ADAQ_data_rasch_na_omit),
  itemnumber = 1:4,
  itemdescrip = c("ADAQ.1", "ADAQ.2", "ADAQ.3", "ADAQ.4"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(ADAQ_data_rasch_na_omit),
  itemnumber = 5:7,
  itemdescrip = c("ADAQ.5", "ADAQ.6", "ADAQ.7"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
#PCA of residuals (assessing unidimensionality)
RIpcmPCA(ADAQ_data_rasch_na_omit)
Eigenvalues Proportion of variance
1.71 24.9%
1.39 21.3%
1.21 20.3%
1.01 12.9%
0.87 11.3%
Code
#item restscores
RIrestscore(ADAQ_data_rasch_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
ADAQ.1 0.56 0.61 0.05 0.414 -0.10 -1.22
ADAQ.2_revised 0.73 0.61 0.12 0.048 * -0.05 -1.17
ADAQ.3 0.65 0.62 0.03 0.701 0.18 -0.95
ADAQ.4 0.54 0.62 0.08 0.414 -0.15 -1.27
ADAQ.5 0.58 0.63 0.05 0.509 -0.01 -1.14
ADAQ.6 0.74 0.64 0.10 0.079 . 0.26 -0.87
ADAQ.7 0.65 0.63 0.02 0.701 -0.14 -1.26
Code
#loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(ADAQ_data_rasch_na_omit)

Code
#residual correlations
resid_cor_ADAQ <- RIgetResidCor(ADAQ_data_rasch_na_omit, iterations = 800, cpu = 10)

RIresidcorr(ADAQ_data_rasch_na_omit, cutoff = resid_cor_ADAQ$p99)
ADAQ.1 ADAQ.2_revised ADAQ.3 ADAQ.4 ADAQ.5 ADAQ.6 ADAQ.7
ADAQ.1
ADAQ.2_revised -0.05
ADAQ.3 -0.26 -0.14
ADAQ.4 -0.23 -0.26 0.1
ADAQ.5 -0.23 0.01 -0.2 -0.08
ADAQ.6 -0.14 -0.08 -0.12 -0.31 -0.25
ADAQ.7 -0.16 -0.11 -0.17 -0.25 -0.36 0.16
Note:
Relative cut-off value is 0.166, which is 0.314 above the average correlation (-0.149).
Correlations above the cut-off are highlighted in red text.
Code
#partial gamma coefficients (assessing local independence)
RIpartgamLD(ADAQ_data_rasch_na_omit) #show only significant results
Item 1 Item 2 Partial gamma SE Lower CI Upper CI Adjusted p-value (BH)
ADAQ.6 ADAQ.7 0.504 0.136 0.238 0.77 0.008
Code
partgam_LD(as.data.frame(ADAQ_data_rasch_na_omit)) #show all output, also non-significant results
            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1          ADAQ.1 ADAQ.2_revised  0.0282 0.1726 0.8703  1.0000     -0.3101
2          ADAQ.1         ADAQ.3 -0.0492 0.1709 0.7735  1.0000     -0.3840
3          ADAQ.1         ADAQ.4 -0.0417 0.1591 0.7934  1.0000     -0.3534
4          ADAQ.1         ADAQ.5 -0.0191 0.1579 0.9037  1.0000     -0.3286
5          ADAQ.1         ADAQ.6  0.0222 0.1783 0.9008  1.0000     -0.3273
6          ADAQ.1         ADAQ.7  0.1077 0.1561 0.4903  1.0000     -0.1983
7  ADAQ.2_revised         ADAQ.3 -0.0377 0.2004 0.8506  1.0000     -0.4304
8  ADAQ.2_revised         ADAQ.4  0.0606 0.1836 0.7414  1.0000     -0.2993
9  ADAQ.2_revised         ADAQ.5  0.3566 0.1599 0.0257  1.0000      0.0433
10 ADAQ.2_revised         ADAQ.6  0.1818 0.1649 0.2703  1.0000     -0.1415
11 ADAQ.2_revised         ADAQ.7  0.1562 0.1735 0.3678  1.0000     -0.1838
12         ADAQ.3         ADAQ.4  0.2436 0.1711 0.1546  1.0000     -0.0918
13         ADAQ.3         ADAQ.5 -0.1408 0.1503 0.3488  1.0000     -0.4355
14         ADAQ.3         ADAQ.6 -0.0145 0.1762 0.9345  1.0000     -0.3599
15         ADAQ.3         ADAQ.7 -0.1220 0.1794 0.4967  1.0000     -0.4736
16         ADAQ.4         ADAQ.5  0.0296 0.1511 0.8447  1.0000     -0.2665
17         ADAQ.4         ADAQ.6 -0.3611 0.1513 0.0170  0.7135     -0.6576
18         ADAQ.4         ADAQ.7 -0.2687 0.1607 0.0945  1.0000     -0.5835
19         ADAQ.5         ADAQ.6 -0.1294 0.1451 0.3724  1.0000     -0.4138
20         ADAQ.5         ADAQ.7 -0.3333 0.1656 0.0441  1.0000     -0.6579
21         ADAQ.6         ADAQ.7  0.5041 0.1356 0.0002  0.0084  **  0.2384
     upper
1   0.3664
2   0.2857
3   0.2701
4   0.2904
5   0.3717
6   0.4137
7   0.3549
8   0.4205
9   0.6699
10  0.5051
11  0.4963
12  0.5790
13  0.1538
14  0.3309
15  0.2297
16  0.3257
17 -0.0646
18  0.0462
19  0.1549
20 -0.0088
21  0.7699

            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1  ADAQ.2_revised         ADAQ.1  0.1241 0.1719 0.4704  1.0000     -0.2129
2          ADAQ.3         ADAQ.1 -0.0850 0.1727 0.6226  1.0000     -0.4234
3          ADAQ.3 ADAQ.2_revised -0.1163 0.1999 0.5608  1.0000     -0.5081
4          ADAQ.4         ADAQ.1 -0.0909 0.1611 0.5726  1.0000     -0.4067
5          ADAQ.4 ADAQ.2_revised -0.1586 0.1792 0.3762  1.0000     -0.5099
6          ADAQ.4         ADAQ.3  0.3691 0.1550 0.0172  0.7242      0.0653
7          ADAQ.5         ADAQ.1  0.0059 0.1527 0.9691  1.0000     -0.2934
8          ADAQ.5 ADAQ.2_revised  0.2500 0.1550 0.1067  1.0000     -0.0537
9          ADAQ.5         ADAQ.3 -0.2000 0.1542 0.1946  1.0000     -0.5022
10         ADAQ.5         ADAQ.4  0.1226 0.1559 0.4318  1.0000     -0.1831
11         ADAQ.6         ADAQ.1  0.1343 0.1662 0.4188  1.0000     -0.1913
12         ADAQ.6 ADAQ.2_revised  0.1966 0.1785 0.2708  1.0000     -0.1533
13         ADAQ.6         ADAQ.3 -0.0484 0.1940 0.8031  1.0000     -0.4287
14         ADAQ.6         ADAQ.4 -0.2667 0.1764 0.1307  1.0000     -0.6125
15         ADAQ.6         ADAQ.5  0.1765 0.1728 0.3072  1.0000     -0.1622
16         ADAQ.7         ADAQ.1 -0.0440 0.1662 0.7911  1.0000     -0.3697
17         ADAQ.7 ADAQ.2_revised -0.0769 0.1810 0.6708  1.0000     -0.4317
18         ADAQ.7         ADAQ.3 -0.1450 0.1775 0.4139  1.0000     -0.4930
19         ADAQ.7         ADAQ.4 -0.1656 0.1712 0.3335  1.0000     -0.5011
20         ADAQ.7         ADAQ.5 -0.2448 0.1613 0.1291  1.0000     -0.5608
21         ADAQ.7         ADAQ.6  0.4437 0.1470 0.0025  0.1070      0.1555
    upper
1  0.4610
2  0.2534
3  0.2755
4  0.2248
5  0.1927
6  0.6729
7  0.3053
8  0.5537
9  0.1022
10 0.4282
11 0.4600
12 0.5464
13 0.3319
14 0.0792
15 0.5152
16 0.2817
17 0.2778
18 0.2029
19 0.1700
20 0.0713
21 0.7319
Code
#targeting
RItargeting(ADAQ_data_rasch_na_omit)

Code
#item hierarchy
itemlabels <- data.frame(itemnr = names(ADAQ_data_rasch_na_omit), item = NA)

RIitemHierarchy(ADAQ_data_rasch_na_omit)

Code
#threshold locations 
RIitemparams(ADAQ_data_rasch_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Threshold 6 Item location
ADAQ.1 -3.04 -2.07 -1.30 1.21 1.38 3.25 -0.1
ADAQ.2_revised -2.62 -1.61 0.56 0.45 2.97 NA -0.05
ADAQ.3 -2.08 -1.94 -1.23 1.30 1.07 3.94 0.18
ADAQ.4 -3.26 -1.66 -1.24 0.94 0.90 3.44 -0.15
ADAQ.5 -2.44 -1.69 -0.39 1.20 0.54 2.69 -0.01
ADAQ.6 -1.48 -0.42 -1.18 1.22 0.55 2.87 0.26
ADAQ.7 -2.02 -2.42 -0.01 0.61 0.48 2.53 -0.14
Note:
Item location is the average of the thresholds for each item.
Code
#analysis of response categories
RIitemCats(ADAQ_data_rasch_na_omit, items = "all", xlims = c(-5, 5), legend = "left")

3.2 Recoding response categories for ADAQ

Based on the analysis of the response categories and threshold locations, scale steps were merged and recoded so that each kept scale step had the highest response probability at some point on the person location scale, and also so that all threshold locations were in the right order for all items. Recoded items were given new names with the abbrevation “_recoded”.

Code
#add recoded ADAQ items 


bulls_eye_merged_superfinal_data <- bulls_eye_merged_superfinal_data %>%
  mutate(
    # Recoding ADAQ.1 (0 + 1)
    ADAQ.1_recoded = recode(
      ADAQ.1,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2,
      `4` = 3,
      `5` = 4,
      `6` = 5
    ),
    # Recoding ADAQ.2 (0 + 1 + 2; 3 + 4)
    ADAQ.2_recoded = recode(
      ADAQ.2,
      `0` = 0, `1` = 0, `2` = 0,
      `3` = 1, `4` = 1,
      `5` = 2,
      `6` = 3
    ),
    # Recoding ADAQ.3 (0 + 1; 3 + 4)
    ADAQ.3_recoded = recode(
      ADAQ.3,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2, `4` = 2,
      `5` = 3,
      `6` = 4
    ),
    
    # Recoding ADAQ.4 (0 + 1)
    ADAQ.4_recoded = recode(
      ADAQ.4,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2,
      `4` = 3,
      `5` = 4,
      `6` = 5
    ),
    
    # Recoding ADAQ.5 (0 + 1; 3 + 4)
    ADAQ.5_recoded = recode(
      ADAQ.5,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2, `4` = 2,
      `5` = 3,
      `6` = 4
    ),
    
    # Recoding ADAQ.6 (#ADAQ.6: 0 + 1 + 2; 3 + 4)
    ADAQ.6_recoded = recode(
      ADAQ.6,
      `0` = 0, `1` = 0, `2` = 0,
      `3` = 1, `4` = 1,
      `5` = 2,
      `6` = 3
    ),
    
    # Recoding ADAQ.7 (#ADAQ.7: 0 + 1, 3 + 4)
    ADAQ.7_recoded = recode(
      ADAQ.7,
      `0` = 0, `1` = 0,
      `2` = 1,
      `3` = 2, `4` = 2,
      `5` = 3,
      `6` = 4
    )
  )

3.3 Re-analysis of ADAQ recoded version

Code
#preparing data frame for ADAQ recoded 
ADAQ_data_rasch_recoded <- bulls_eye_merged_superfinal_data[, c(
  "ADAQ.1_recoded",
  "ADAQ.2_recoded",
  "ADAQ.3_recoded",
  "ADAQ.4_recoded",
  "ADAQ.5_recoded",
  "ADAQ.6_recoded",
  "ADAQ.7_recoded"
)]

#investigating missing values
RImissing(ADAQ_data_rasch_recoded)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete ADAQ 
ADAQ_data_rasch_recoded_na_omit <- ADAQ_data_rasch_recoded %>%
  select(ADAQ.1_recoded, ADAQ.2_recoded, ADAQ.3_recoded, ADAQ.4_recoded, ADAQ.5_recoded, ADAQ.6_recoded, ADAQ.7_recoded) %>% 
  na.omit()

RImissing(ADAQ_data_rasch_recoded_na_omit)
[1] "No missing data."
Code
#re-analysis of response categories using ADAQ recoded
RIitemCats(ADAQ_data_rasch_recoded_na_omit, items = "all", xlims = c(-8, 8), legend = "left")

Code
#tile plot (response category frequency for the recoded items) 
RItileplot(ADAQ_data_rasch_recoded_na_omit)

Code
#re-analysis reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI
RItif(na.omit(ADAQ_data_rasch_recoded_na_omit))

Code
RItif(na.omit(ADAQ_data_rasch_recoded_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(ADAQ_data_rasch_recoded_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr df pvalue sig
overall 26  27 0.52
Code
#re-analysis of conditional item fit (assessing unidimensionaliity)
simfit1_ADAQ_recoded <- RIgetfit(ADAQ_data_rasch_recoded_na_omit, iterations = 100, cpu = 10) 

RIitemfit(ADAQ_data_rasch_recoded_na_omit, simfit1_ADAQ_recoded, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
ADAQ.1_recoded 1.416 [0.759, 1.241] 1.478 [0.74, 1.26] 0.175 0.218 -0.67
ADAQ.2_recoded 0.716 [0.817, 1.215] 0.708 [0.811, 1.272] 0.101 0.103 -0.37
ADAQ.3_recoded 0.771 [0.775, 1.305] 0.741 [0.775, 1.387] 0.004 0.034 -0.57
ADAQ.4_recoded 1.355 [0.713, 1.303] 1.449 [0.683, 1.332] 0.052 0.117 -0.69
ADAQ.5_recoded 1.087 [0.736, 1.237] 1.107 [0.707, 1.343] no misfit no misfit -0.71
ADAQ.6_recoded 0.726 [0.758, 1.216] 0.693 [0.755, 1.359] 0.032 0.062 0.14
ADAQ.7_recoded 0.962 [0.754, 1.314] 0.958 [0.751, 1.347] no misfit no misfit -0.95
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 97 simulated datasets.
Code
#re-analysis of conditional item charachteristic curves
ICCplot(
  as.data.frame(ADAQ_data_rasch_recoded_na_omit),
  itemnumber = 1:4,
  itemdescrip = c(
    "ADAQ.1_recoded",
    "ADAQ.2_recoded",
    "ADAQ.3_recoded",
    "ADAQ.4_recoded"
  ),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(ADAQ_data_rasch_recoded_na_omit),
  itemnumber = 5:7,
  itemdescrip = c("ADAQ.5_recoded", "ADAQ.6_recoded", "ADAQ.7_recoded"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
#re-analyzing of PCA of residuals (assessing unidimensionality)
RIpcmPCA(ADAQ_data_rasch_recoded_na_omit)
Eigenvalues Proportion of variance
1.66 26%
1.40 23.5%
1.20 18.7%
1.10 13.3%
0.87 9.7%
Code
#re-analysis of item restscores
RIrestscore(ADAQ_data_rasch_recoded_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
ADAQ.1_recoded 0.54 0.67 0.13 0.032 * -0.08 -0.67
ADAQ.2_recoded 0.78 0.63 0.15 0.003 ** 0.21 -0.37
ADAQ.3_recoded 0.75 0.64 0.11 0.033 * 0.02 -0.57
ADAQ.4_recoded 0.55 0.67 0.12 0.090 . -0.11 -0.69
ADAQ.5_recoded 0.62 0.65 0.03 0.623 -0.12 -0.71
ADAQ.6_recoded 0.76 0.64 0.12 0.032 * 0.72 0.14
ADAQ.7_recoded 0.68 0.65 0.03 0.623 -0.37 -0.95
Code
#re-analysis of loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(ADAQ_data_rasch_recoded_na_omit)

Code
#re-analysis of residual correlations
resid_cor_ADAQ_recoded <- RIgetResidCor(ADAQ_data_rasch_recoded_na_omit,
                                        iterations = 800,
                                        cpu = 10)

RIresidcorr(ADAQ_data_rasch_recoded_na_omit, cutoff = resid_cor_ADAQ_recoded$p99)
ADAQ.1_recoded ADAQ.2_recoded ADAQ.3_recoded ADAQ.4_recoded ADAQ.5_recoded ADAQ.6_recoded ADAQ.7_recoded
ADAQ.1_recoded
ADAQ.2_recoded -0.16
ADAQ.3_recoded -0.29 0.08
ADAQ.4_recoded -0.28 -0.27 -0.01
ADAQ.5_recoded -0.27 -0.03 -0.11 -0.1
ADAQ.6_recoded -0.2 0.03 0 -0.34 -0.12
ADAQ.7_recoded -0.17 -0.04 -0.14 -0.28 -0.31 0.2
Note:
Relative cut-off value is 0.173, which is 0.307 above the average correlation (-0.134).
Correlations above the cut-off are highlighted in red text.
Code
#re-analysis of partial gamma coefficients (assessing local independence)
RIpartgamLD(ADAQ_data_rasch_recoded_na_omit) #shows only significant results
Item 1 Item 2 Partial gamma SE Lower CI Upper CI Adjusted p-value (BH)
ADAQ.6_recoded ADAQ.7_recoded 0.54 0.138 0.269 0.811 0.004
Code
partgam_LD(as.data.frame(ADAQ_data_rasch_recoded_na_omit)) #shows all output, also non-significant results
            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1  ADAQ.1_recoded ADAQ.2_recoded -0.0345 0.1860 0.8529  1.0000     -0.3990
2  ADAQ.1_recoded ADAQ.3_recoded -0.2952 0.1828 0.1063  1.0000     -0.6535
3  ADAQ.1_recoded ADAQ.4_recoded -0.0816 0.1491 0.5840  1.0000     -0.3739
4  ADAQ.1_recoded ADAQ.5_recoded  0.0000 0.1461 1.0000  1.0000     -0.2863
5  ADAQ.1_recoded ADAQ.6_recoded -0.0541 0.1993 0.7862  1.0000     -0.4447
6  ADAQ.1_recoded ADAQ.7_recoded  0.1011 0.1483 0.4954  1.0000     -0.1896
7  ADAQ.2_recoded ADAQ.3_recoded  0.3256 0.2527 0.1976  1.0000     -0.1696
8  ADAQ.2_recoded ADAQ.4_recoded  0.0794 0.2173 0.7150  1.0000     -0.3466
9  ADAQ.2_recoded ADAQ.5_recoded  0.3162 0.1934 0.1021  1.0000     -0.0629
10 ADAQ.2_recoded ADAQ.6_recoded  0.2778 0.2247 0.2163  1.0000     -0.1626
11 ADAQ.2_recoded ADAQ.7_recoded  0.2182 0.2095 0.2976  1.0000     -0.1923
12 ADAQ.3_recoded ADAQ.4_recoded  0.1533 0.2054 0.4555  1.0000     -0.2493
13 ADAQ.3_recoded ADAQ.5_recoded -0.0515 0.2015 0.7981  1.0000     -0.4465
14 ADAQ.3_recoded ADAQ.6_recoded -0.0309 0.2559 0.9038  1.0000     -0.5324
15 ADAQ.3_recoded ADAQ.7_recoded -0.0957 0.2121 0.6520  1.0000     -0.5113
16 ADAQ.4_recoded ADAQ.5_recoded  0.0307 0.1596 0.8476  1.0000     -0.2821
17 ADAQ.4_recoded ADAQ.6_recoded -0.3617 0.1684 0.0317  1.0000     -0.6917
18 ADAQ.4_recoded ADAQ.7_recoded -0.2067 0.1697 0.2233  1.0000     -0.5394
19 ADAQ.5_recoded ADAQ.6_recoded -0.1240 0.2048 0.5450  1.0000     -0.5254
20 ADAQ.5_recoded ADAQ.7_recoded -0.3780 0.1647 0.0217  0.9111     -0.7008
21 ADAQ.6_recoded ADAQ.7_recoded  0.5397 0.1383 0.0001  0.0040  **  0.2686
     upper
1   0.3300
2   0.0630
3   0.2106
4   0.2863
5   0.3366
6   0.3918
7   0.8208
8   0.5053
9   0.6954
10  0.7181
11  0.6287
12  0.5558
13  0.3434
14  0.4705
15  0.3200
16  0.3435
17 -0.0317
18  0.1260
19  0.2775
20 -0.0553
21  0.8107

            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1  ADAQ.2_recoded ADAQ.1_recoded  0.2319 0.1989 0.2437  1.0000     -0.1580
2  ADAQ.3_recoded ADAQ.1_recoded -0.3097 0.1877 0.0989  1.0000     -0.6776
3  ADAQ.3_recoded ADAQ.2_recoded  0.3253 0.2403 0.1758  1.0000     -0.1456
4  ADAQ.4_recoded ADAQ.1_recoded -0.1337 0.1815 0.4614  1.0000     -0.4894
5  ADAQ.4_recoded ADAQ.2_recoded -0.3200 0.1677 0.0564  1.0000     -0.6487
6  ADAQ.4_recoded ADAQ.3_recoded  0.1875 0.1931 0.3316  1.0000     -0.1910
7  ADAQ.5_recoded ADAQ.1_recoded -0.0508 0.1726 0.7682  1.0000     -0.3890
8  ADAQ.5_recoded ADAQ.2_recoded  0.1579 0.2008 0.4316  1.0000     -0.2356
9  ADAQ.5_recoded ADAQ.3_recoded -0.1376 0.2001 0.4916  1.0000     -0.5297
10 ADAQ.5_recoded ADAQ.4_recoded  0.1088 0.1509 0.4710  1.0000     -0.1870
11 ADAQ.6_recoded ADAQ.1_recoded  0.1971 0.1933 0.3080  1.0000     -0.1818
12 ADAQ.6_recoded ADAQ.2_recoded  0.2039 0.2296 0.3745  1.0000     -0.2461
13 ADAQ.6_recoded ADAQ.3_recoded -0.0787 0.2740 0.7741  1.0000     -0.6158
14 ADAQ.6_recoded ADAQ.4_recoded -0.2222 0.1873 0.2356  1.0000     -0.5894
15 ADAQ.6_recoded ADAQ.5_recoded -0.0392 0.2119 0.8532  1.0000     -0.4545
16 ADAQ.7_recoded ADAQ.1_recoded  0.1638 0.1602 0.3064  1.0000     -0.1501
17 ADAQ.7_recoded ADAQ.2_recoded  0.0492 0.2003 0.8060  1.0000     -0.3433
18 ADAQ.7_recoded ADAQ.3_recoded -0.1667 0.2121 0.4319  1.0000     -0.5823
19 ADAQ.7_recoded ADAQ.4_recoded -0.0296 0.1703 0.8621  1.0000     -0.3634
20 ADAQ.7_recoded ADAQ.5_recoded -0.2893 0.1828 0.1135  1.0000     -0.6475
21 ADAQ.7_recoded ADAQ.6_recoded  0.4468 0.1765 0.0114  0.4773      0.1008
    upper
1  0.6218
2  0.0581
3  0.7962
4  0.2220
5  0.0087
6  0.5660
7  0.2873
8  0.5514
9  0.2545
10 0.4046
11 0.5760
12 0.6539
13 0.4585
14 0.1450
15 0.3760
16 0.4778
17 0.4417
18 0.2490
19 0.3042
20 0.0689
21 0.7928
Code
#re-analysis of targeting
RItargeting(ADAQ_data_rasch_recoded_na_omit)

Code
#re-analysis of item hierarchy
itemlabels <- data.frame(itemnr = names(ADAQ_data_rasch_recoded_na_omit), item = NA)

RIitemHierarchy(ADAQ_data_rasch_recoded_na_omit)

Code
#threshold locations 
RIitemparams(ADAQ_data_rasch_recoded_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Item location
ADAQ.1_recoded -3.39 -2.05 0.82 1.17 3.05 -0.08
ADAQ.2_recoded -2.72 0.64 2.72 NA NA 0.21
ADAQ.3_recoded -2.89 -2.17 1.50 3.65 NA 0.02
ADAQ.4_recoded -2.92 -2.01 0.52 0.65 3.23 -0.11
ADAQ.5_recoded -2.61 -1.28 0.98 2.43 NA -0.12
ADAQ.6_recoded -1.41 0.97 2.61 NA NA 0.72
ADAQ.7_recoded -3.31 -1.10 0.66 2.28 NA -0.37
Note:
Item location is the average of the thresholds for each item.

Removing item ADAQ.6_recoded as it does not work well and is potentially redundant. Items ADAQ.1_recoded, ADAQ.2_recoded, and ADAQ.4_recoded does not perform that well neither but is kept for now in the next analysis step. ADAQ.6_recoded is the item that performs worst.

3.4 Re-analysis of ADAQ recoded and trimmed version

The ADAQ.6_recoded item is dropped as it did not perform well and is potentially redundant. A re-analysis is conducted with the recoded items but without ADAQ.6_recoded, called the trimmed version.

Code
#preparing data frame for ADAQ recoded and trimmed version

ADAQ_data_rasch_recoded_trimmed <- bulls_eye_merged_superfinal_data[, c(
  "ADAQ.1_recoded",
  "ADAQ.2_recoded",
  "ADAQ.3_recoded",
  "ADAQ.4_recoded",
  "ADAQ.5_recoded",
  "ADAQ.7_recoded"
)]

#investigating missing values
RImissing(ADAQ_data_rasch_recoded_trimmed)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete ADAQ
ADAQ_data_rasch_recoded_trimmed_na_omit <- ADAQ_data_rasch_recoded_trimmed %>%
  select(ADAQ.1_recoded, ADAQ.2_recoded, ADAQ.3_recoded, ADAQ.4_recoded, ADAQ.5_recoded, ADAQ.7_recoded) %>% 
  na.omit()

RImissing(ADAQ_data_rasch_recoded_trimmed_na_omit)
[1] "No missing data."
Code
#re-analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI
RItif(na.omit(ADAQ_data_rasch_recoded_trimmed_na_omit))

Code
RItif(na.omit(ADAQ_data_rasch_recoded_trimmed_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(ADAQ_data_rasch_recoded_trimmed_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr df pvalue sig
overall 20  24 0.7
Code
#re-analysis of conditional item fit (assessing unidimensionaliity)
simfit1_ADAQ_recoded_trimmed <- RIgetfit(ADAQ_data_rasch_recoded_trimmed_na_omit, iterations = 100, cpu = 10) 

RIitemfit(ADAQ_data_rasch_recoded_trimmed_na_omit, simfit1_ADAQ_recoded_trimmed, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
ADAQ.1_recoded 1.337 [0.785, 1.305] 1.365 [0.777, 1.266] 0.032 0.099 -0.60
ADAQ.2_recoded 0.709 [0.806, 1.257] 0.7 [0.789, 1.26] 0.097 0.089 -0.35
ADAQ.3_recoded 0.746 [0.766, 1.223] 0.72 [0.702, 1.23] 0.02 no misfit -0.49
ADAQ.4_recoded 1.22 [0.708, 1.264] 1.282 [0.712, 1.278] no misfit 0.004 -0.62
ADAQ.5_recoded 1.028 [0.777, 1.228] 1.037 [0.773, 1.281] no misfit no misfit -0.64
ADAQ.7_recoded 0.992 [0.83, 1.28] 0.999 [0.822, 1.298] no misfit no misfit -0.88
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 96 simulated datasets.
Code
#re-analysis of conditional item charachteristic curves
ICCplot(
  as.data.frame(ADAQ_data_rasch_recoded_trimmed_na_omit),
  itemnumber = 1:4,
  itemdescrip = c(
    "ADAQ.1_recoded",
    "ADAQ.2_recoded",
    "ADAQ.3_recoded",
    "ADAQ.4_recoded"
  )
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(ADAQ_data_rasch_recoded_trimmed_na_omit),
  itemnumber = 5:6,
  itemdescrip = c("ADAQ.5_recoded", "ADAQ.7_recoded")
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
#re-analyzing PCA of residuals (assessing unidimensionality)
RIpcmPCA(ADAQ_data_rasch_recoded_trimmed_na_omit)
Eigenvalues Proportion of variance
1.49 28.8%
1.41 24.2%
1.21 20.9%
1.05 15.2%
0.83 10.8%
Code
#re-analysis of item restscores
RIrestscore(ADAQ_data_rasch_recoded_trimmed_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
ADAQ.1_recoded 0.53 0.64 0.11 0.066 . 0.01 -0.60
ADAQ.2_recoded 0.77 0.61 0.16 0.001 ** 0.26 -0.35
ADAQ.3_recoded 0.75 0.61 0.14 0.017 * 0.12 -0.49
ADAQ.4_recoded 0.56 0.65 0.09 0.298 -0.01 -0.62
ADAQ.5_recoded 0.62 0.63 0.01 0.901 -0.03 -0.64
ADAQ.7_recoded 0.64 0.62 0.02 0.901 -0.27 -0.88
Code
#re-analysis of loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(ADAQ_data_rasch_recoded_trimmed_na_omit)

Code
#re-analysis of residual correlations
resid_cor_ADAQ_recoded_trimmed <- RIgetResidCor(ADAQ_data_rasch_recoded_trimmed_na_omit,
                                                iterations = 800,
                                                cpu = 10)

RIresidcorr(ADAQ_data_rasch_recoded_trimmed_na_omit, cutoff = resid_cor_ADAQ_recoded_trimmed$p99)
ADAQ.1_recoded ADAQ.2_recoded ADAQ.3_recoded ADAQ.4_recoded ADAQ.5_recoded ADAQ.7_recoded
ADAQ.1_recoded
ADAQ.2_recoded -0.16
ADAQ.3_recoded -0.3 0.11
ADAQ.4_recoded -0.35 -0.3 -0.02
ADAQ.5_recoded -0.3 -0.02 -0.1 -0.14
ADAQ.7_recoded -0.15 0 -0.09 -0.28 -0.27
Note:
Relative cut-off value is 0.144, which is 0.303 above the average correlation (-0.159).
Correlations above the cut-off are highlighted in red text.
Code
#re-analysis of partial gamma coefficients (assessing local independence)
RIpartgamLD(ADAQ_data_rasch_recoded_trimmed_na_omit) #shows only significant results
[1] "No statistically significant local dependency found."
Code
partgam_LD(as.data.frame(ADAQ_data_rasch_recoded_trimmed_na_omit)) #shows all output, also non-significant results
            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1  ADAQ.1_recoded ADAQ.2_recoded -0.0122 0.1796 0.9459  1.0000     -0.3641
2  ADAQ.1_recoded ADAQ.3_recoded -0.2441 0.1829 0.1820  1.0000     -0.6025
3  ADAQ.1_recoded ADAQ.4_recoded -0.1250 0.1561 0.4232  1.0000     -0.4309
4  ADAQ.1_recoded ADAQ.5_recoded -0.0200 0.1526 0.8957  1.0000     -0.3190
5  ADAQ.1_recoded ADAQ.7_recoded  0.1429 0.1507 0.3432  1.0000     -0.1525
6  ADAQ.2_recoded ADAQ.3_recoded  0.3167 0.1967 0.1075  1.0000     -0.0689
7  ADAQ.2_recoded ADAQ.4_recoded -0.1176 0.1816 0.5171  1.0000     -0.4736
8  ADAQ.2_recoded ADAQ.5_recoded  0.3951 0.1691 0.0195  0.5837      0.0637
9  ADAQ.2_recoded ADAQ.7_recoded  0.3333 0.1725 0.0533  1.0000     -0.0048
10 ADAQ.3_recoded ADAQ.4_recoded  0.3000 0.1714 0.0801  1.0000     -0.0360
11 ADAQ.3_recoded ADAQ.5_recoded -0.0161 0.2069 0.9379  1.0000     -0.4216
12 ADAQ.3_recoded ADAQ.7_recoded  0.0946 0.2005 0.6371  1.0000     -0.2984
13 ADAQ.4_recoded ADAQ.5_recoded -0.1810 0.1731 0.2959  1.0000     -0.5202
14 ADAQ.4_recoded ADAQ.7_recoded -0.2444 0.1684 0.1466  1.0000     -0.5745
15 ADAQ.5_recoded ADAQ.7_recoded -0.2077 0.1775 0.2420  1.0000     -0.5557
    upper
1  0.3398
2  0.1144
3  0.1809
4  0.2790
5  0.4383
6  0.7022
7  0.2383
8  0.7264
9  0.6715
10 0.6360
11 0.3893
12 0.4876
13 0.1583
14 0.0856
15 0.1402

            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1  ADAQ.2_recoded ADAQ.1_recoded  0.3168 0.1727 0.0666  1.0000     -0.0217
2  ADAQ.3_recoded ADAQ.1_recoded -0.1719 0.2070 0.4064  1.0000     -0.5777
3  ADAQ.3_recoded ADAQ.2_recoded  0.3274 0.1995 0.1007  1.0000     -0.0635
4  ADAQ.4_recoded ADAQ.1_recoded -0.0977 0.1660 0.5563  1.0000     -0.4231
5  ADAQ.4_recoded ADAQ.2_recoded -0.4356 0.1512 0.0040  0.1186     -0.7318
6  ADAQ.4_recoded ADAQ.3_recoded  0.0350 0.1938 0.8568  1.0000     -0.3449
7  ADAQ.5_recoded ADAQ.1_recoded -0.0802 0.1496 0.5919  1.0000     -0.3733
8  ADAQ.5_recoded ADAQ.2_recoded  0.1975 0.1829 0.2804  1.0000     -0.1611
9  ADAQ.5_recoded ADAQ.3_recoded -0.1603 0.1908 0.4007  1.0000     -0.5342
10 ADAQ.5_recoded ADAQ.4_recoded -0.0239 0.1506 0.8737  1.0000     -0.3190
11 ADAQ.7_recoded ADAQ.1_recoded  0.1152 0.1574 0.4642  1.0000     -0.1933
12 ADAQ.7_recoded ADAQ.2_recoded  0.1818 0.1874 0.3321  1.0000     -0.1856
13 ADAQ.7_recoded ADAQ.3_recoded  0.0483 0.2191 0.8256  1.0000     -0.3812
14 ADAQ.7_recoded ADAQ.4_recoded -0.0101 0.1712 0.9529  1.0000     -0.3456
15 ADAQ.7_recoded ADAQ.5_recoded -0.0909 0.1986 0.6471  1.0000     -0.4801
     upper
1   0.6552
2   0.2339
3   0.7184
4   0.2277
5  -0.1393
6   0.4148
7   0.2129
8   0.5560
9   0.2136
10  0.2712
11  0.4237
12  0.5492
13  0.4778
14  0.3254
15  0.2983
Code
#re-analysis of targeting
RItargeting(ADAQ_data_rasch_recoded_trimmed_na_omit)

Code
#re-analysis of item hierarchy
itemlabels <- data.frame(itemnr = names(ADAQ_data_rasch_recoded_trimmed_na_omit),
                         item = NA)

RIitemHierarchy(ADAQ_data_rasch_recoded_trimmed_na_omit)

Code
#threshold locations 
RIitemparams(ADAQ_data_rasch_recoded_trimmed_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Item location
ADAQ.1_recoded -3.01 -1.86 0.84 1.12 2.95 0.01
ADAQ.2_recoded -2.49 0.65 2.64 NA NA 0.26
ADAQ.3_recoded -2.51 -2.02 1.47 3.53 NA 0.12
ADAQ.4_recoded -2.55 -1.82 0.54 0.62 3.13 -0.01
ADAQ.5_recoded -2.27 -1.17 0.96 2.35 NA -0.03
ADAQ.7_recoded -2.96 -0.99 0.65 2.21 NA -0.27
Note:
Item location is the average of the thresholds for each item.
Code
RIitemCats(ADAQ_data_rasch_recoded_trimmed_na_omit, items = "all", xlims = c(-8, 8), legend = "left")

Items ADAQ.1_recoded, ADAQ.2_recoded, and ADAQ.3_recoded still does not perform optimal (see conditional item fit and item restscores) but is kept in the scale for the validity analysis.

3.5 ADAQ and BEA validity analysis

Convergent validity between ADAQ and BEA is investigated by correlational analysis. As both scales has undergone Rasch analysis, participants theta estimates (computed with weighted likelihood estimation - WLE) will be correlated, and since thetas are interval data Pearson’s correlation is applied.

3.6 Data preparation for validity analyses

The following data frame will be used for all validity analyses.

Code
#complete cases data frame for validity analyses

#removing missing values in relevant items for complete cases analysis
#only BE items and one item from the other scales in the bulls eye study
#is needed (in this case the SubjektivPrestation item) 
#this is because there are no missing values for the other scales (beside some BE items) 
#after the intervention study (Hockeystudien) participants have been removed
#therefore this data frame will work for all validity analyses 

items_validity_analyses <- c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded",
  "SubjektivPrestation"
)

data_complete_cases <- bulls_eye_merged_superfinal_data[complete.cases(bulls_eye_merged_superfinal_data[, items_validity_analyses]), ]

#156 participants to begin with
#34 participants from the intervention study (Hockeystudien) were removed (did not complete any scale beside Bulls Eye)
#8 participants (not from the intervention study) with missing values in Bulls Eye items were removed
#results in a complete cases data frame with 114 participants 

complete_cases_dims <- dim_names(data_complete_cases, "Användarnamn") 
count <- dplyr::count
count(complete_cases_dims)
# A tibble: 1 × 1
      n
  <int>
1   114
Code
#item data frames
BE_validity_analyses <- data_complete_cases[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

ADAQ_validity_analyses <- data_complete_cases[, c(
  "ADAQ.1_recoded",
  "ADAQ.2_recoded",
  "ADAQ.3_recoded",
  "ADAQ.4_recoded",
  "ADAQ.5_recoded",
  "ADAQ.7_recoded"
)]

#extracting thetas
BE_thetas <- RIestThetas(BE_validity_analyses)

ADAQ_thetas <- RIestThetas(ADAQ_validity_analyses)

#adding thetas to the data frame
data_complete_cases$BE_thetas_WLE <- BE_thetas$WLE
data_complete_cases$ADAQ_thetas_WLE <- ADAQ_thetas$WLE
Code
correlation_BE_ADAQ <- cor_test(
  data_complete_cases,
  x = "BE_thetas_WLE",
  y = "ADAQ_thetas_WLE",
  method = "pearson"
)

print(correlation_BE_ADAQ)
Parameter1    |      Parameter2 |    r |       95% CI | t(112) |         p
--------------------------------------------------------------------------
BE_thetas_WLE | ADAQ_thetas_WLE | 0.33 | [0.16, 0.49] |   3.71 | < .001***

Observations: 114
Code
plot(correlation_BE_ADAQ)

3.6.1 Sensitivity analysis ADAQ and BEA, time latency

Some participants did not complete BEA and ADAQ during the same day as intended. Therefore, a sensitivity analysis is conducted.

Code
#time latency variables have previously been created
#using dataframe from correlation calculation to account for removed NA
mean(data_complete_cases$time_latency_survey1_bulls_eye1_days,
     na.rm = TRUE)
[1] 0.877088
Code
sd(data_complete_cases$time_latency_survey1_bulls_eye1_days,
   na.rm = TRUE)
[1] 2.824687
Code
range(data_complete_cases$time_latency_survey1_bulls_eye1_days,
      na.rm = TRUE)
[1]  0.003171296 21.032523148
Code
#count participants with a time latency of more than 1 day 
count_greater_than_1 <- sum(data_complete_cases$time_latency_survey1_bulls_eye1_days > 1,
                            na.rm = TRUE)

#count participants with a time latency of more than 3 days 
count_greater_than_3 <- sum(data_complete_cases$time_latency_survey1_bulls_eye1_days > 3,
                            na.rm = TRUE)

print(paste("Number of participants with a time latency of more than 1 day:", count_greater_than_1))
[1] "Number of participants with a time latency of more than 1 day: 15"
Code
print(paste("Number of participants with a time latency of more than 3 days:", count_greater_than_3))
[1] "Number of participants with a time latency of more than 3 days: 10"
Code
#creating data frame with only participants that completed Bulls Eye and the other scales during the same day
#since 15 participants did not, 99 participants will remain in the sensitivity analysis
data_complete_cases_less_than_1day <- subset(data_complete_cases, time_latency_survey1_bulls_eye1_days < 1)
Code
#sensitivity analysis ADAQ and BE
BE_validity_analyses_less1day <- data_complete_cases_less_than_1day[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

ADAQ_validity_analyses_less1day <- data_complete_cases_less_than_1day[, c(
  "ADAQ.1_recoded",
  "ADAQ.2_recoded",
  "ADAQ.3_recoded",
  "ADAQ.4_recoded",
  "ADAQ.5_recoded",
  "ADAQ.7_recoded"
)]

#extracting thetas
BE_thetas_less1day <- RIestThetas(BE_validity_analyses_less1day)

ADAQ_thetas_less1day <- RIestThetas(ADAQ_validity_analyses_less1day)

#adding thetas to the data frame
data_complete_cases_less_than_1day$BE_thetas_less1day_WLE <- BE_thetas_less1day$WLE
data_complete_cases_less_than_1day$ADAQ_thetas_less1day_WLE <- ADAQ_thetas_less1day$WLE

#analysis only with participants completing Bulls Eye and ADAQ the same day
correlation_BE_ADAQ_less1day <- cor_test(
  data_complete_cases_less_than_1day,
  x = "BE_thetas_less1day_WLE",
  y = "ADAQ_thetas_less1day_WLE",
  method = "pearson"
)

print(correlation_BE_ADAQ_less1day)
Parameter1             |               Parameter2 |    r |       95% CI | t(97) |         p
-------------------------------------------------------------------------------------------
BE_thetas_less1day_WLE | ADAQ_thetas_less1day_WLE | 0.34 | [0.15, 0.50] |  3.58 | < .001***

Observations: 99
Code
#we can also check if there is a statistically significant difference 
#by adding same or different day as a predictor 
data_complete_cases_sensitivity_analysis <- data_complete_cases %>% 
  mutate(group = factor(ifelse(time_latency_survey1_bulls_eye1_days < 1,"same_day","different_day")))

m0_adaq <- lm(ADAQ_thetas_WLE ~ BE_thetas_WLE, data = data_complete_cases_sensitivity_analysis)
m1_adaq <- lm(ADAQ_thetas_WLE ~ BE_thetas_WLE + group, data = data_complete_cases_sensitivity_analysis)

anova(m0_adaq, m1_adaq)
Analysis of Variance Table

Model 1: ADAQ_thetas_WLE ~ BE_thetas_WLE
Model 2: ADAQ_thetas_WLE ~ BE_thetas_WLE + group
  Res.Df    RSS Df Sum of Sq      F Pr(>F)
1    112 273.64                           
2    111 270.25  1    3.3951 1.3945 0.2402

Adding whether BE and ADAQ were completed during the same day or not as a predictor did not significantly alter the results.

4 Satisfaction With Life Scale (SWLS): Rasch analysis and BEA validity investigation

The SWLS instrument is used as part of the validity investigation of BEA. First, SWLS will undergo Rasch analysis. Then a validity analysis of correlating BEA and SWLS thetas will be undertaken. Note that SWLS was not completed by the participants from the intervention study (Hockeystudien), therefore more missing values than in the BEA Rasch analysis are to be expected.

4.1 Rasch analysis of SWLS

Code
#preparing data frame with SWLS items 
SWLS_data_rasch <- bulls_eye_merged_superfinal_data[, c("SWLS.1", "SWLS.2", "SWLS.3", "SWLS.4", "SWLS.5")]


#investigate missing values
RImissing(SWLS_data_rasch)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete SWLS
SWLS_data_rasch_na_omit <- SWLS_data_rasch %>%
  select(SWLS.1, SWLS.2, SWLS.3, SWLS.4, SWLS.5) %>% 
  na.omit()

RImissing(SWLS_data_rasch_na_omit)
[1] "No missing data."
Code
#tile plot (response category frequency for the items)
RItileplot(SWLS_data_rasch_na_omit)

SWLS.2 has no responses in scale step one and will therefore be merged with scale step zero. SWLS.3 and SWLS.4 has no responses in scale step zero and will therefore be merged with scale step one.

Code
#merging scale step zero and one in SWLS.2, SWLS.3, and SWLS.4 (creating SWLS.2_revised, SWLS.3_revised, and SWLS.4_revised)
bulls_eye_merged_superfinal_data <- bulls_eye_merged_superfinal_data %>%
  mutate(
    # Recoding SWLS.2 (0+1)
    SWLS.2_revised = recode(SWLS.2,
                            `0` = 0, `1` = 0, 
                            `2` = 1,
                            `3` = 2, 
                            `4` = 3,
                            `5` = 4,
                            `6` = 5),
    # Recoding SWLS.3 (0+1)
    SWLS.3_revised = recode(SWLS.3,
                            `0` = 0, `1` = 0, 
                            `2` = 1,
                            `3` = 2, 
                            `4` = 3,
                            `5` = 4,
                            `6` = 5),
    # Recoding SWLS.4 (0+1)
    SWLS.4_revised = recode(SWLS.4,
                            `0` = 0, `1` = 0, 
                            `2` = 1,
                            `3` = 2, 
                            `4` = 3,
                            `5` = 4,
                            `6` = 5))


#preparing data frame with SWLS items including SWLS.2_revised, SWLS.3_revised, and SWLS.4_revised
SWLS_data_rasch <- bulls_eye_merged_superfinal_data[, c("SWLS.1", "SWLS.2_revised", "SWLS.3_revised", "SWLS.4_revised", "SWLS.5")]


#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete SWLS
SWLS_data_rasch_na_omit <- SWLS_data_rasch %>%
  select(SWLS.1, SWLS.2_revised, SWLS.3_revised, SWLS.4_revised, SWLS.5) %>% 
  na.omit()

RImissing(SWLS_data_rasch_na_omit)
[1] "No missing data."
Code
RItileplot(SWLS_data_rasch_na_omit)

Code
#analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI 
RItif(na.omit(SWLS_data_rasch_na_omit))

Code
RItif(na.omit(SWLS_data_rasch_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(SWLS_data_rasch_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue sig
overall 15.8 26 0.94
Code
#conditional item fit (assessing unidimensionaliity)
simfit1_SWLS <- RIgetfit(SWLS_data_rasch_na_omit, iterations = 100, cpu = 10) 

RIitemfit(SWLS_data_rasch_na_omit, simfit1_SWLS, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
SWLS.1 1.096 [0.8, 1.222] 1.213 [0.843, 1.203] no misfit 0.01 -0.48
SWLS.2_revised 1.128 [0.79, 1.203] 1.188 [0.763, 1.28] no misfit no misfit -0.87
SWLS.3_revised 0.814 [0.779, 1.285] 0.714 [0.713, 1.39] no misfit no misfit -0.86
SWLS.4_revised 0.775 [0.787, 1.196] 0.839 [0.759, 1.266] 0.012 no misfit -0.93
SWLS.5 1.2 [0.735, 1.279] 1.162 [0.734, 1.285] no misfit no misfit -0.32
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 80 simulated datasets.
Code
#conditional item charachteristic curves
ICCplot(
  as.data.frame(SWLS_data_rasch_na_omit),
  itemnumber = 1:4,
  itemdescrip = c("SWLS.1", "SWLS.2", "SWLS.3", "SWLS.4"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(SWLS_data_rasch_na_omit),
  itemnumber = 5:5,
  itemdescrip = c("SWLS.5"),
  method = "cut"
)

Code
#PCA of residuals (assessing unidimensionality)
RIpcmPCA(SWLS_data_rasch_na_omit)
Eigenvalues Proportion of variance
1.49 30.4%
1.34 29.5%
1.27 26%
0.89 14%
0.01 0.1%
Code
#item restscore
RIrestscore(SWLS_data_rasch_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
SWLS.1 0.44 0.50 0.06 0.700 0.19 -0.48
SWLS.2_revised 0.47 0.50 0.03 0.749 -0.20 -0.87
SWLS.3_revised 0.64 0.50 0.14 0.083 . -0.19 -0.86
SWLS.4_revised 0.57 0.49 0.08 0.615 -0.26 -0.93
SWLS.5 0.50 0.52 0.02 0.749 0.35 -0.32
Code
#loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(SWLS_data_rasch_na_omit)

Code
#residual correlations
resid_cor_SWLS <- RIgetResidCor(SWLS_data_rasch_na_omit,
                                iterations = 800,
                                cpu = 10)

RIresidcorr(SWLS_data_rasch_na_omit, cutoff = resid_cor_SWLS$p99)
SWLS.1 SWLS.2_revised SWLS.3_revised SWLS.4_revised SWLS.5
SWLS.1
SWLS.2_revised -0.24
SWLS.3_revised -0.12 -0.23
SWLS.4_revised -0.17 -0.13 0.12
SWLS.5 -0.33 -0.26 -0.28 -0.28
Note:
Relative cut-off value is 0.141, which is 0.333 above the average correlation (-0.192).
Correlations above the cut-off are highlighted in red text.
Code
#partial gamma coefficients (assessing local independence)
RIpartgamLD(SWLS_data_rasch_na_omit) #shows only significant results
Item 1 Item 2 Partial gamma SE Lower CI Upper CI Adjusted p-value (BH)
SWLS.3_revised SWLS.4_revised 0.47 0.149 0.178 0.762 0.032
Code
partgam_LD(as.data.frame(SWLS_data_rasch_na_omit)) #shows all output, also non-significant results
            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1          SWLS.1 SWLS.2_revised -0.0373 0.1498 0.8032  1.0000     -0.3310
2          SWLS.1 SWLS.3_revised -0.2050 0.1631 0.2088  1.0000     -0.5247
3          SWLS.1 SWLS.4_revised -0.1682 0.1411 0.2333  1.0000     -0.4449
4          SWLS.1         SWLS.5 -0.1139 0.1370 0.4057  1.0000     -0.3825
5  SWLS.2_revised SWLS.3_revised -0.1416 0.1678 0.3988  1.0000     -0.4705
6  SWLS.2_revised SWLS.4_revised -0.0448 0.1521 0.7684  1.0000     -0.3429
7  SWLS.2_revised         SWLS.5  0.0615 0.1384 0.6568  1.0000     -0.2097
8  SWLS.3_revised SWLS.4_revised  0.4703 0.1489 0.0016  0.0318   *  0.1784
9  SWLS.3_revised         SWLS.5  0.1808 0.1491 0.2252  1.0000     -0.1114
10 SWLS.4_revised         SWLS.5 -0.0563 0.1498 0.7068  1.0000     -0.3499
    upper
1  0.2563
2  0.1147
3  0.1084
4  0.1546
5  0.1873
6  0.2533
7  0.3327
8  0.7622
9  0.4730
10 0.2372

            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1  SWLS.2_revised         SWLS.1 -0.0120 0.1343 0.9285       1     -0.2753
2  SWLS.3_revised         SWLS.1  0.0588 0.1791 0.7426       1     -0.2922
3  SWLS.3_revised SWLS.2_revised -0.0652 0.1765 0.7118       1     -0.4112
4  SWLS.4_revised         SWLS.1 -0.1897 0.1686 0.2606       1     -0.5201
5  SWLS.4_revised SWLS.2_revised  0.0093 0.1615 0.9543       1     -0.3074
6  SWLS.4_revised SWLS.3_revised  0.2511 0.1606 0.1180       1     -0.0637
7          SWLS.5         SWLS.1  0.0574 0.1271 0.6518       1     -0.1918
8          SWLS.5 SWLS.2_revised  0.0871 0.1479 0.5559       1     -0.2028
9          SWLS.5 SWLS.3_revised  0.1250 0.1575 0.4275       1     -0.1838
10         SWLS.5 SWLS.4_revised -0.0968 0.1484 0.5144       1     -0.3877
    upper
1  0.2512
2  0.4099
3  0.2807
4  0.1407
5  0.3259
6  0.5659
7  0.3066
8  0.3771
9  0.4338
10 0.1941
Code
#targeting
RItargeting(SWLS_data_rasch_na_omit)

Code
#item hierarchy
itemlabels <- data.frame(itemnr = names(SWLS_data_rasch_na_omit), item = NA)

RIitemHierarchy(SWLS_data_rasch_na_omit)

Code
#threshold locations 
RIitemparams(SWLS_data_rasch_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Threshold 6 Item location
SWLS.1 -1.81 -1.13 0.40 -0.59 1.22 3.05 0.19
SWLS.2_revised -2.27 -0.01 -0.90 0.61 1.57 NA -0.2
SWLS.3_revised -1.79 0.47 -1.46 0.33 1.52 NA -0.19
SWLS.4_revised -1.81 -0.71 -0.89 0.54 1.55 NA -0.26
SWLS.5 -0.23 -0.80 0.58 -0.35 0.83 2.06 0.35
Note:
Item location is the average of the thresholds for each item.
Code
#analysis of response categories
RIitemCats(SWLS_data_rasch_na_omit, items = "all", xlims = c(-5, 5), legend = "left")

4.2 Recoding response categories for SWLS

Based on the analysis of the response categories and threshold locations, scale steps were merged and recoded so that each kept scale step had the highest response probability at some point on the person location scale, and also so that all threshold locations were in the right order for all items. Recoded items were given new names with the abbrevation “_recoded”.

Code
#add recoded SWLS items 

bulls_eye_merged_superfinal_data <- bulls_eye_merged_superfinal_data %>%
  mutate(
    # Recoding SWLS.1 (0 + 1; 2 + 3)
    SWLS.1_recoded = recode(SWLS.1,
                            `0` = 0, `1` = 0,
                            `2` = 1, `3` = 1,
                            `4` = 2,
                            `5` = 3,
                            `6` = 4),
    # Recoding SWLS.2 (0+1+2; 3+4)
    SWLS.2_recoded = recode(SWLS.2,
                            `0` = 0, `1` = 0, `2` = 0,
                            `3` = 1, `4` = 1,
                            `5` = 2,
                            `6` = 3),
    # Recoding SWLS.3 (0+1+2; 3+4)
    SWLS.3_recoded = recode(SWLS.3,
                            `0` = 0, `1` = 0, `2` = 0,
                            `3` = 1, `4` = 1,
                            `5` = 2,
                            `6` = 3),
    # Recoding SWLS.4 (0+1+2+3)
    SWLS.4_recoded = recode(SWLS.4,
                            `0` = 0, `1` = 0, `2` = 0, `3` = 0,
                            `4` = 1,
                            `5` = 2,
                            `6` = 3),
    # Recoding SWLS.5 (0+1, 2+3)
    SWLS.5_recoded = recode(SWLS.5,
                            `0` = 0, `1` = 0,
                            `2` = 1, `3` = 1,
                            `4` = 2,
                            `5` = 3,
                            `6` = 4))

4.3 Re-analysis of SWLS recoded version

Code
#preparing data frame for SWLS recoded
SWLS_data_rasch_recoded <- bulls_eye_merged_superfinal_data[, c(
  "SWLS.1_recoded",
  "SWLS.2_recoded",
  "SWLS.3_recoded",
  "SWLS.4_recoded",
  "SWLS.5_recoded"
)]

#investigate missing values 
RImissing(SWLS_data_rasch_recoded)

Code
#removing missing values
SWLS_data_rasch_recoded_na_omit <- SWLS_data_rasch_recoded %>%
  select(SWLS.1_recoded, SWLS.2_recoded, SWLS.3_recoded, SWLS.4_recoded, SWLS.5_recoded) %>% 
  na.omit()

RImissing(SWLS_data_rasch_recoded_na_omit)
[1] "No missing data."
Code
#re-analysis of response categories using SWLS recoded
RIitemCats(SWLS_data_rasch_recoded_na_omit, items = "all", xlims = c(-8, 8), legend = "left")

Code
#tile plot (response category frequency for the recoded items) 
RItileplot(SWLS_data_rasch_recoded_na_omit)

Code
#re-analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI 
RItif(na.omit(SWLS_data_rasch_recoded_na_omit))

Code
RItif(na.omit(SWLS_data_rasch_recoded_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(SWLS_data_rasch_recoded_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue sig
overall 14.4 16 0.57
Code
#re-analysis of conditional item fit (assessing unidimensionaliity)
simfit1_SWLS_recoded <- RIgetfit(SWLS_data_rasch_recoded_na_omit, iterations = 100, cpu = 10) 

RIitemfit(SWLS_data_rasch_recoded_na_omit, simfit1_SWLS_recoded, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
SWLS.1_recoded 1.146 [0.731, 1.229] 1.154 [0.718, 1.231] no misfit no misfit 0.14
SWLS.2_recoded 1.09 [0.801, 1.193] 1.173 [0.761, 1.367] no misfit no misfit -0.29
SWLS.3_recoded 0.769 [0.809, 1.245] 0.725 [0.742, 1.424] 0.04 0.017 -0.40
SWLS.4_recoded 0.851 [0.749, 1.311] 0.853 [0.699, 1.298] no misfit no misfit -0.16
SWLS.5_recoded 1.122 [0.773, 1.297] 1.186 [0.739, 1.386] no misfit no misfit 0.00
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 100 simulated datasets.
Code
#re-analysis of conditional item charachteristic curves
ICCplot(
  as.data.frame(SWLS_data_rasch_recoded_na_omit),
  itemnumber = 1:4,
  itemdescrip = c(
    "SWLS.1_recoded",
    "SWLS.2_recoded",
    "SWLS.3_recoded",
    "SWLS.4_recoded"
  ),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(SWLS_data_rasch_recoded_na_omit),
  itemnumber = 5,
  itemdescrip = c("SWLS.5_recoded"),
  method = "cut"
)

Code
#re-analysis of PCA of residuals (assessing unidimensionality)
RIpcmPCA(SWLS_data_rasch_recoded_na_omit)
Eigenvalues Proportion of variance
1.47 30.1%
1.33 29.2%
1.28 25.7%
0.91 14.8%
0.01 0.1%
Code
#re-analysis of item restscores
RIrestscore(SWLS_data_rasch_recoded_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
SWLS.1_recoded 0.44 0.52 0.08 0.441 0.25 0.14
SWLS.2_recoded 0.48 0.52 0.04 0.775 -0.18 -0.29
SWLS.3_recoded 0.66 0.52 0.14 0.095 . -0.28 -0.40
SWLS.4_recoded 0.60 0.52 0.08 0.441 -0.04 -0.16
SWLS.5_recoded 0.52 0.54 0.02 0.775 0.12 0.00
Code
#re-analyzing loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(SWLS_data_rasch_recoded_na_omit)

Code
#re-analysis of residual correlations
resid_cor_SWLS_recoded <- RIgetResidCor(SWLS_data_rasch_recoded_na_omit,
                                        iterations = 800,
                                        cpu = 10)

RIresidcorr(SWLS_data_rasch_recoded_na_omit, cutoff = resid_cor_SWLS_recoded$p99)
SWLS.1_recoded SWLS.2_recoded SWLS.3_recoded SWLS.4_recoded SWLS.5_recoded
SWLS.1_recoded
SWLS.2_recoded -0.25
SWLS.3_recoded -0.17 -0.26
SWLS.4_recoded -0.2 -0.13 0.09
SWLS.5_recoded -0.3 -0.21 -0.18 -0.36
Note:
Relative cut-off value is 0.091, which is 0.289 above the average correlation (-0.198).
Correlations above the cut-off are highlighted in red text.
Code
#re-analysis of partial gamma coefficients (assessing local independence)
RIpartgamLD(SWLS_data_rasch_recoded_na_omit) #shows only significant results
[1] "No statistically significant local dependency found."
Code
partgam_LD(as.data.frame(SWLS_data_rasch_recoded_na_omit)) #shows all output, also non-significant results
            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1  SWLS.1_recoded SWLS.2_recoded -0.1184 0.1573 0.4518       1     -0.4267
2  SWLS.1_recoded SWLS.3_recoded -0.0855 0.1608 0.5949       1     -0.4006
3  SWLS.1_recoded SWLS.4_recoded -0.1408 0.1350 0.2970       1     -0.4054
4  SWLS.1_recoded SWLS.5_recoded -0.1090 0.1435 0.4474       1     -0.3903
5  SWLS.2_recoded SWLS.3_recoded -0.2996 0.1735 0.0842       1     -0.6396
6  SWLS.2_recoded SWLS.4_recoded -0.0244 0.1595 0.8785       1     -0.3371
7  SWLS.2_recoded SWLS.5_recoded  0.0389 0.1666 0.8155       1     -0.2877
8  SWLS.3_recoded SWLS.4_recoded  0.3278 0.1745 0.0603       1     -0.0142
9  SWLS.3_recoded SWLS.5_recoded  0.2555 0.1764 0.1475       1     -0.0902
10 SWLS.4_recoded SWLS.5_recoded -0.1357 0.1568 0.3867       1     -0.4430
    upper
1  0.1900
2  0.2296
3  0.1238
4  0.1722
5  0.0405
6  0.2883
7  0.3654
8  0.6698
9  0.6012
10 0.1716

            Item1          Item2   gamma     se pvalue padj.BH sig   lower
1  SWLS.2_recoded SWLS.1_recoded -0.0817 0.1531 0.5935       1     -0.3818
2  SWLS.3_recoded SWLS.1_recoded  0.1312 0.1626 0.4196       1     -0.1874
3  SWLS.3_recoded SWLS.2_recoded -0.1429 0.1962 0.4665       1     -0.5274
4  SWLS.4_recoded SWLS.1_recoded  0.0678 0.1708 0.6914       1     -0.2669
5  SWLS.4_recoded SWLS.2_recoded  0.0311 0.1788 0.8618       1     -0.3193
6  SWLS.4_recoded SWLS.3_recoded  0.1508 0.2009 0.4529       1     -0.2429
7  SWLS.5_recoded SWLS.1_recoded -0.0347 0.1392 0.8029       1     -0.3077
8  SWLS.5_recoded SWLS.2_recoded  0.1048 0.1572 0.5047       1     -0.2032
9  SWLS.5_recoded SWLS.3_recoded  0.1111 0.1775 0.5313       1     -0.2367
10 SWLS.5_recoded SWLS.4_recoded -0.1698 0.1484 0.2524       1     -0.4606
    upper
1  0.2183
2  0.4499
3  0.2417
4  0.4025
5  0.3815
6  0.5445
7  0.2382
8  0.4129
9  0.4590
10 0.1210
Code
#re-analysis of targeting
RItargeting(SWLS_data_rasch_recoded_na_omit)

Code
#re-analysis of item hierarchy
itemlabels <- data.frame(itemnr = names(SWLS_data_rasch_recoded_na_omit), item = NA)

RIitemHierarchy(SWLS_data_rasch_recoded_na_omit)

Code
#threshold locations 
RIitemparams(SWLS_data_rasch_recoded_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Item location
SWLS.1_recoded -2.10 -0.53 0.64 3.02 0.25
SWLS.2_recoded -1.84 0.15 1.16 NA -0.18
SWLS.3_recoded -1.69 -0.22 1.07 NA -0.28
SWLS.4_recoded -1.14 -0.13 1.15 NA -0.04
SWLS.5_recoded -1.26 -0.23 0.17 1.80 0.12
Note:
Item location is the average of the thresholds for each item.

SWLS.3 was a candidate for removal but after adjusting response categories it performed slightly better (at least in the gamma coefficient analysis), so it is kept in the scale for the validity analysis.

4.4 SWLS and BEA validity analysis

Validity between SWLS and BEA is investigated by correlational analysis. As both scales has undergone Rasch analysis, participants theta estimates (computed with weighted likelihood estimation - WLE) will be correlated, and since thetas are interval data Pearson’s correlation is applied.

Code
#complete cases data frame for validity analyses

#removing missing values in relevant items for complete cases analysis
#only BE items and one item from the other scales in the bulls eye study
#is needed (in this case the SubjektivPrestation item) 
#this is because there are no missing values for the other scales (beside some BE items) 
#after the intervention study (Hockeystudien) participants have been removed
#therefore this data frame will work for all validity analyses 

items_validity_analyses <- c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded",
  "SubjektivPrestation"
)

data_complete_cases <- bulls_eye_merged_superfinal_data[complete.cases(bulls_eye_merged_superfinal_data[, items_validity_analyses]), ]

#156 participants to begin with
#34 participants from the intervention study (Hockeystudien) were removed (did not complete any scale beside Bulls Eye)
#8 participants (not from the intervention study) with missing values in Bulls Eye items were removed
#results in a complete cases data frame with 114 participants 

complete_cases_dims <- dim_names(data_complete_cases, "Användarnamn") 
count <- dplyr::count
count(complete_cases_dims)
# A tibble: 1 × 1
      n
  <int>
1   114
Code
#item data frames
BE_validity_analyses <- data_complete_cases[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

SWLS_validity_analyses <- data_complete_cases[, c(
  "SWLS.1_recoded",
  "SWLS.2_recoded",
  "SWLS.3_recoded",
  "SWLS.4_recoded",
  "SWLS.5_recoded"
)]

#extracting thetas
BE_thetas <- RIestThetas(BE_validity_analyses)

SWLS_thetas <- RIestThetas(SWLS_validity_analyses)

#adding thetas to the data frame
data_complete_cases$BE_thetas_WLE <- BE_thetas$WLE 
data_complete_cases$SWLS_thetas_WLE <- SWLS_thetas$WLE
Code
correlation_BE_SWLS <- cor_test(
  data_complete_cases,
  x = "BE_thetas_WLE",
  y = "SWLS_thetas_WLE",
  method = "pearson"
)

print(correlation_BE_SWLS)
Parameter1    |      Parameter2 |    r |       95% CI | t(112) |         p
--------------------------------------------------------------------------
BE_thetas_WLE | SWLS_thetas_WLE | 0.57 | [0.43, 0.68] |   7.35 | < .001***

Observations: 114
Code
plot(correlation_BE_SWLS)

4.4.1 Sensitivity analysis SWLS and BEA, time latency

Some participants did not complete BEA and SWLS during the same day as intended. Therefore, a sensitivity analysis is conducted.

Code
#creating data frame with only participants that completed Bulls Eye and the other scales during the same day
#since 15 participants did not, 99 participants will remain in the sensitivity analysis
data_complete_cases_less_than_1day <- subset(data_complete_cases, time_latency_survey1_bulls_eye1_days < 1)
Code
#data have previously been prepared for the validity sensitivity analysis
#looking only at participants who completed Bulls Eye and the other scales during the same day (N = 99)
#BE items and BE thetas for these particpants have also been prepared previously

#sensitivity analysis SWLS and BE
BE_validity_analyses_less1day <- data_complete_cases_less_than_1day[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

SWLS_validity_analyses_less1day <- data_complete_cases_less_than_1day[, c(
  "SWLS.1_recoded",
  "SWLS.2_recoded",
  "SWLS.3_recoded",
  "SWLS.4_recoded",
  "SWLS.5_recoded"
)]

#extracting thetas
BE_thetas_less1day <- RIestThetas(BE_validity_analyses_less1day)

SWLS_thetas_less1day <- RIestThetas(SWLS_validity_analyses_less1day)

#adding thetas to the data frame
data_complete_cases_less_than_1day$BE_thetas_less1day_WLE <- BE_thetas_less1day$WLE
data_complete_cases_less_than_1day$SWLS_thetas_less1day_WLE <- SWLS_thetas_less1day$WLE

#analysis only with participants completing Bulls Eye and ADAQ the same day
correlation_BE_SWLS_less1day <- cor_test(
  data_complete_cases_less_than_1day,
  x = "BE_thetas_less1day_WLE",
  y = "SWLS_thetas_less1day_WLE",
  method = "pearson"
)

print(correlation_BE_SWLS_less1day)
Parameter1             |               Parameter2 |    r |       95% CI | t(97) |         p
-------------------------------------------------------------------------------------------
BE_thetas_less1day_WLE | SWLS_thetas_less1day_WLE | 0.56 | [0.41, 0.68] |  6.63 | < .001***

Observations: 99
Code
#we can also check if there is a statistically significant difference 
#by adding same or different day as a predictor 
data_complete_cases_sensitivity_analysis <- data_complete_cases %>% 
  mutate(group = factor(ifelse(time_latency_survey1_bulls_eye1_days < 1,"same_day","different_day")))

m0_swls <- lm(SWLS_thetas_WLE ~ BE_thetas_WLE, data = data_complete_cases_sensitivity_analysis)
m1_swls <- lm(SWLS_thetas_WLE ~ BE_thetas_WLE + group, data = data_complete_cases_sensitivity_analysis)

anova(m0_swls, m1_swls)
Analysis of Variance Table

Model 1: SWLS_thetas_WLE ~ BE_thetas_WLE
Model 2: SWLS_thetas_WLE ~ BE_thetas_WLE + group
  Res.Df    RSS Df Sum of Sq      F Pr(>F)
1    112 137.03                           
2    111 137.00  1  0.023626 0.0191 0.8902

Adding whether BE and SWLS were completed during the same day or not as a predictor did not significantly alter the results.

5 Valuing Questionnaire (VQ): Rasch analysis and BEA validity investigation

The Valuing Questionnaire (VQ) instrument is used as part of the validity investigation of BEA. The original VQ study (Smout et al., 2014) states that the scale has two dimensions: (A) progress in valued living (items 3, 4, 5, 7, 9), and (B) obstruction to valued living (items 1, 2, 6, 8, 10). VQ_progress and VQ_obstruction are therefore analyzed separate with Rasch analysis. Validity analyses by correlations of thetas (BEA and VQ_progress, and BEA and VQ_obstruction respectively) will be undertaken. Note that VQ was not completed by the participants from the intervention study (Hockeystudien), therefore more missing values than in the BEA Rasch analysis are to be expected.

#VQ_progess: Rasch analysis and BEA validity investigation

Code
#preparing data frame with VQ_progress items
VQ_progress_data_rasch <- bulls_eye_merged_superfinal_data[, c("VQ.3", "VQ.4", "VQ.5", "VQ.7", "VQ.9")]

#investigate missing values
RImissing(VQ_progress_data_rasch)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete VQ
VQ_progress_data_rasch_na_omit <- VQ_progress_data_rasch %>%
  select(VQ.3,
         VQ.4,
         VQ.5,
         VQ.7,
         VQ.9) %>%
  na.omit()

RImissing(VQ_progress_data_rasch_na_omit)
[1] "No missing data."
Code
#tile plot (response category frequency for the items)
RItileplot(VQ_progress_data_rasch_na_omit)

Code
#analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI 
RItif(na.omit(VQ_progress_data_rasch_na_omit))

Code
RItif(na.omit(VQ_progress_data_rasch_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(VQ_progress_data_rasch_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue  sig 
overall 61.4 29 0.00042  ***
Code
#conditional item fit (assessing unidimensionaliity)
simfit1_VQ_progress <- RIgetfit(VQ_progress_data_rasch_na_omit,
                       iterations = 100,
                       cpu = 10) 

RIitemfit(VQ_progress_data_rasch_na_omit, simfit1_VQ_progress, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
VQ.3 2.058 [0.754, 1.253] 2.519 [0.73, 1.263] 0.805 1.256 -0.26
VQ.4 0.818 [0.766, 1.276] 0.757 [0.742, 1.286] no misfit no misfit -0.72
VQ.5 0.731 [0.814, 1.254] 0.764 [0.817, 1.211] 0.083 0.053 -0.37
VQ.7 0.698 [0.81, 1.237] 0.72 [0.789, 1.258] 0.112 0.069 -0.75
VQ.9 0.777 [0.777, 1.167] 0.812 [0.714, 1.184] no misfit no misfit -0.71
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 63 simulated datasets.
Code
#conditional item charachteristic curves
ICCplot(
  as.data.frame(VQ_progress_data_rasch_na_omit),
  itemnumber = 1:4,
  itemdescrip = c("VQ.3", "VQ.4", "VQ.5", "VQ.7"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(VQ_progress_data_rasch_na_omit),
  itemnumber = 5,
  itemdescrip = c("VQ.9"),
  method = "cut"
)

Code
#PCA of residuals (assessing unidimensionality)
RIpcmPCA(VQ_progress_data_rasch_na_omit)
Eigenvalues Proportion of variance
1.87 55%
1.26 18.2%
1.08 15.6%
0.78 11%
0.01 0.2%
Code
#item restscore
RIrestscore(VQ_progress_data_rasch_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
VQ.3 0.12 0.45 0.33 0.000 *** 0.30 -0.26
VQ.4 0.56 0.44 0.12 0.075 . -0.15 -0.72
VQ.5 0.58 0.44 0.14 0.033 * 0.19 -0.37
VQ.7 0.61 0.42 0.19 0.001 *** -0.19 -0.75
VQ.9 0.56 0.43 0.13 0.070 . -0.15 -0.71
Code
#loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(VQ_progress_data_rasch_na_omit)

Code
#residual correlations
resid_cor_VQ_progress <- RIgetResidCor(VQ_progress_data_rasch_na_omit,
                                       iterations = 800,
                                       cpu = 10)

RIresidcorr(VQ_progress_data_rasch_na_omit, cutoff = resid_cor_VQ_progress$p99)
VQ.3 VQ.4 VQ.5 VQ.7 VQ.9
VQ.3
VQ.4 -0.6
VQ.5 -0.35 -0.07
VQ.7 -0.4 0.22 0.1
VQ.9 -0.41 0.14 -0.02 -0.08
Note:
Relative cut-off value is 0.159, which is 0.306 above the average correlation (-0.146).
Correlations above the cut-off are highlighted in red text.
Code
#partial gamma coefficients (assessing local independence)
RIpartgamLD(VQ_progress_data_rasch_na_omit) #shows only significant results
Item 1 Item 2 Partial gamma SE Lower CI Upper CI Adjusted p-value (BH)
VQ.7 VQ.4 0.391 0.128 0.141 0.642 0.044
VQ.7 VQ.5 0.368 0.108 0.156 0.580 0.013
Code
partgam_LD(as.data.frame(VQ_progress_data_rasch_na_omit)) #shows all output, also non-significant results
   Item1 Item2   gamma     se pvalue padj.BH  sig   lower   upper
1   VQ.3  VQ.4 -0.5596 0.1072 0.0000  0.0000  *** -0.7697 -0.3494
2   VQ.3  VQ.5 -0.2258 0.1125 0.0448  0.8965      -0.4464 -0.0052
3   VQ.3  VQ.7 -0.4400 0.1284 0.0006  0.0122    * -0.6916 -0.1884
4   VQ.3  VQ.9 -0.3835 0.1252 0.0022  0.0437    * -0.6288 -0.1382
5   VQ.4  VQ.5  0.1811 0.1387 0.1915  1.0000      -0.0907  0.4529
6   VQ.4  VQ.7  0.4236 0.1443 0.0033  0.0668    .  0.1407  0.7065
7   VQ.4  VQ.9  0.3306 0.1311 0.0116  0.2329       0.0738  0.5875
8   VQ.5  VQ.7  0.3054 0.1415 0.0309  0.6188       0.0280  0.5829
9   VQ.5  VQ.9  0.1368 0.1366 0.3169  1.0000      -0.1310  0.4045
10  VQ.7  VQ.9  0.1598 0.1649 0.3323  1.0000      -0.1633  0.4829

   Item1 Item2   gamma     se pvalue padj.BH sig   lower   upper
1   VQ.4  VQ.3 -0.3189 0.1342 0.0174  0.3490     -0.5819 -0.0560
2   VQ.5  VQ.3  0.2364 0.1216 0.0519  1.0000     -0.0020  0.4747
3   VQ.5  VQ.4  0.1970 0.1374 0.1516  1.0000     -0.0722  0.4662
4   VQ.7  VQ.3  0.0000 0.1503 1.0000  1.0000     -0.2947  0.2947
5   VQ.7  VQ.4  0.3913 0.1278 0.0022  0.0439   *  0.1409  0.6417
6   VQ.7  VQ.5  0.3684 0.1081 0.0007  0.0132   *  0.1565  0.5804
7   VQ.9  VQ.3  0.0294 0.1459 0.8402  1.0000     -0.2565  0.3153
8   VQ.9  VQ.4  0.3252 0.1296 0.0121  0.2426      0.0711  0.5793
9   VQ.9  VQ.5  0.0315 0.1443 0.8272  1.0000     -0.2513  0.3143
10  VQ.9  VQ.7  0.0171 0.1635 0.9167  1.0000     -0.3033  0.3375
Code
#targeting
RItargeting(VQ_progress_data_rasch_na_omit)

Code
#item hierarchy
itemlabels <- data.frame(itemnr = names(VQ_progress_data_rasch_na_omit), item = NA)

RIitemHierarchy(VQ_progress_data_rasch_na_omit)

Code
#threshold locations 
RIitemparams(VQ_progress_data_rasch_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Threshold 6 Item location
VQ.3 0.23 -0.18 -0.58 -0.03 1.24 1.12 0.3
VQ.4 -1.16 -1.32 0.17 0.06 0.71 0.62 -0.15
VQ.5 -0.64 -0.44 -0.48 0.48 0.81 1.42 0.19
VQ.7 -1.70 -1.10 -0.50 -0.18 0.86 1.48 -0.19
VQ.9 -1.18 -1.08 -0.15 -0.01 0.70 0.84 -0.15
Note:
Item location is the average of the thresholds for each item.
Code
#analysis of response categories
RIitemCats(VQ_progress_data_rasch_na_omit, items = "all", xlims = c(-5, 5), legend = "left")

5.1 Recoding response categories for VQ_progress

Based on the analysis of the response categories and threshold locations, scale steps were merged and recoded so that each kept scale step had the highest response probability at some point on the person location scale, and also so that all threshold locations were in the right order for all items. Recoded items were given new names with the abbrevation “_recoded”.

Code
#add recoded VQ_progress items

bulls_eye_merged_superfinal_data <- bulls_eye_merged_superfinal_data %>%
  mutate(
    # Recoding VQ.3 (0+1; 2+3)
    VQ.3_recoded = recode(VQ.3,
                            `0` = 0, `1` = 0, 
                            `2` = 1, `3` = 1, 
                            `4` = 2,
                            `5` = 3,
                            `6` = 4),
    # Recoding VQ.4 (0+1; 3+4; 5+6)
    VQ.4_recoded = recode(VQ.4,
                            `0` = 0, `1` = 0, 
                            `2` = 1, 
                            `3` = 2, `4` = 2,
                            `5` = 3, `6` = 3),
    # Recoding VQ.5 (0+1+2)
    VQ.5_recoded = recode(VQ.5,
                            `0` = 0, `1` = 0,`2` = 0,
                            `3` = 1,
                            `4` = 2,
                            `5` = 3,
                            `6` = 4),
    # Recoding VQ.7 (0+1)
    VQ.7_recoded = recode(VQ.7,
                           `0` = 0, `1` = 0, 
                           `2` = 1,
                           `3` = 2, 
                           `4` = 3,
                           `5` = 4,
                           `6` = 5),
    # Recoding VQ.9 (0+1)
    VQ.9_recoded = recode(VQ.9,
                           `0` = 0, `1` = 0, 
                           `2` = 1, 
                           `3` = 2, 
                           `4` = 3,
                           `5` = 4, 
                           `6` = 5))

5.2 Re-analysis of VQ_progress recoded version

Code
#preparing items for VQ_progress recoded
VQ_progress_data_rasch_recoded <- bulls_eye_merged_superfinal_data[, c(
  "VQ.3_recoded",
  "VQ.4_recoded",
  "VQ.5_recoded",
  "VQ.7_recoded",
  "VQ.9_recoded"
)]

#investigate missing values
RImissing(VQ_progress_data_rasch_recoded)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete VQ
VQ_progress_data_rasch_recoded_na_omit <- VQ_progress_data_rasch_recoded %>%
  select(
    VQ.3_recoded,
    VQ.4_recoded,
    VQ.5_recoded,
    VQ.7_recoded,
    VQ.9_recoded
  ) %>%
  na.omit()

RImissing(VQ_progress_data_rasch_recoded_na_omit)
[1] "No missing data."
Code
#re-analysis of response categories using VQ_progress recoded
RIitemCats(
  VQ_progress_data_rasch_recoded_na_omit,
  items = "all",
  xlims = c(-5, 5),
  legend = "left"
)

Code
#tile plot (response category frequency for the items)
RItileplot(VQ_progress_data_rasch_recoded_na_omit)

Code
#re-analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI 
RItif(na.omit(VQ_progress_data_rasch_recoded_na_omit))

Code
RItif(na.omit(VQ_progress_data_rasch_recoded_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(VQ_progress_data_rasch_recoded_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue sig
overall 40.8 20 0.004   **
Code
#re-analysis of conditional item fit (assessing unidimensionaliity)
simfit1_VQ_progress_recoded <- RIgetfit(VQ_progress_data_rasch_recoded_na_omit,
                               iterations = 100,
                               cpu = 10)

RIitemfit(VQ_progress_data_rasch_recoded_na_omit, simfit1_VQ_progress_recoded, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
VQ.3_recoded 1.821 [0.792, 1.333] 1.898 [0.758, 1.343] 0.488 0.555 -0.02
VQ.4_recoded 0.721 [0.824, 1.245] 0.657 [0.839, 1.294] 0.103 0.182 -1.00
VQ.5_recoded 0.804 [0.773, 1.262] 0.804 [0.734, 1.258] no misfit no misfit 0.15
VQ.7_recoded 0.762 [0.779, 1.291] 0.768 [0.773, 1.308] 0.017 0.005 -0.47
VQ.9_recoded 0.947 [0.741, 1.308] 0.944 [0.701, 1.413] no misfit no misfit -0.50
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 100 simulated datasets.
Code
#re-analysis of conditional item charachteristic curves
ICCplot(
  as.data.frame(VQ_progress_data_rasch_recoded_na_omit),
  itemnumber = 1:4,
  itemdescrip = c(
    "VQ.3_recoded",
    "VQ.4_recoded",
    "VQ.5_recoded",
    "VQ.7_recoded"
  ),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(VQ_progress_data_rasch_recoded_na_omit),
  itemnumber = 5,
  itemdescrip = c("VQ.9_recoded"),
  method = "cut"
)

Code
#re-analysis of PCA of residuals (assessing unidimensionality)
RIpcmPCA(VQ_progress_data_rasch_recoded_na_omit)
Eigenvalues Proportion of variance
1.74 48.3%
1.37 22.6%
1.11 17.5%
0.76 11.2%
0.01 0.3%
Code
#re-analysis of item restscores
RIrestscore(VQ_progress_data_rasch_recoded_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
VQ.3_recoded 0.13 0.48 0.35 0.000 *** 0.32 -0.02
VQ.4_recoded 0.66 0.43 0.23 0.000 *** -0.65 -1.00
VQ.5_recoded 0.61 0.48 0.13 0.043 * 0.50 0.15
VQ.7_recoded 0.61 0.48 0.13 0.023 * -0.12 -0.47
VQ.9_recoded 0.55 0.49 0.06 0.404 -0.16 -0.50
Code
#re-analysis of loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(VQ_progress_data_rasch_recoded_na_omit)

Code
#re-analysis of residual correlations
resid_cor_VQ_progress_recoded <- RIgetResidCor(VQ_progress_data_rasch_recoded_na_omit,
                                               iterations = 800,
                                               cpu = 10)

RIresidcorr(VQ_progress_data_rasch_recoded_na_omit, cutoff = resid_cor_VQ_progress_recoded$p99)
VQ.3_recoded VQ.4_recoded VQ.5_recoded VQ.7_recoded VQ.9_recoded
VQ.3_recoded
VQ.4_recoded -0.46
VQ.5_recoded -0.3 -0.06
VQ.7_recoded -0.39 0.2 0.04
VQ.9_recoded -0.41 0.13 -0.19 -0.2
Note:
Relative cut-off value is 0.15, which is 0.314 above the average correlation (-0.164).
Correlations above the cut-off are highlighted in red text.
Code
#re-analysis of partial gamma coefficients (assessing local independence)
RIpartgamLD(VQ_progress_data_rasch_recoded_na_omit) #shows only significant results
Item 1 Item 2 Partial gamma SE Lower CI Upper CI Adjusted p-value (BH)
VQ.4_recoded VQ.7_recoded 0.566 0.126 0.320 0.812 0.000
VQ.4_recoded VQ.9_recoded 0.431 0.130 0.176 0.686 0.018
Code
partgam_LD(as.data.frame(VQ_progress_data_rasch_recoded_na_omit)) #shows all output, also non-significant results
          Item1        Item2   gamma     se pvalue padj.BH  sig   lower   upper
1  VQ.3_recoded VQ.4_recoded -0.4865 0.1267 0.0001  0.0025   ** -0.7348 -0.2382
2  VQ.3_recoded VQ.5_recoded -0.2375 0.1233 0.0540  1.0000      -0.4791  0.0041
3  VQ.3_recoded VQ.7_recoded -0.3704 0.1367 0.0067  0.1347      -0.6383 -0.1025
4  VQ.3_recoded VQ.9_recoded -0.2781 0.1264 0.0278  0.5568      -0.5259 -0.0303
5  VQ.4_recoded VQ.5_recoded  0.2165 0.1774 0.2224  1.0000      -0.1313  0.5643
6  VQ.4_recoded VQ.7_recoded  0.5660 0.1256 0.0000  0.0001  ***  0.3198  0.8123
7  VQ.4_recoded VQ.9_recoded  0.4309 0.1299 0.0009  0.0183    *  0.1762  0.6856
8  VQ.5_recoded VQ.7_recoded  0.2285 0.1519 0.1327  1.0000      -0.0693  0.5263
9  VQ.5_recoded VQ.9_recoded  0.0186 0.1430 0.8966  1.0000      -0.2617  0.2989
10 VQ.7_recoded VQ.9_recoded  0.1073 0.1498 0.4740  1.0000      -0.1864  0.4009

          Item1        Item2   gamma     se pvalue padj.BH sig   lower  upper
1  VQ.4_recoded VQ.3_recoded -0.2000 0.1486 0.1783  1.0000     -0.4912 0.0912
2  VQ.5_recoded VQ.3_recoded  0.2000 0.1286 0.1198  1.0000     -0.0520 0.4520
3  VQ.5_recoded VQ.4_recoded  0.0211 0.1989 0.9157  1.0000     -0.3689 0.4110
4  VQ.7_recoded VQ.3_recoded -0.0246 0.1553 0.8742  1.0000     -0.3290 0.2798
5  VQ.7_recoded VQ.4_recoded  0.3679 0.1463 0.0119  0.2384      0.0811 0.6547
6  VQ.7_recoded VQ.5_recoded  0.3498 0.1199 0.0035  0.0705   .  0.1148 0.5848
7  VQ.9_recoded VQ.3_recoded -0.0312 0.1572 0.8424  1.0000     -0.3393 0.2768
8  VQ.9_recoded VQ.4_recoded  0.3237 0.1368 0.0179  0.3588      0.0556 0.5917
9  VQ.9_recoded VQ.5_recoded -0.0758 0.1571 0.6296  1.0000     -0.3836 0.2321
10 VQ.9_recoded VQ.7_recoded  0.0108 0.1554 0.9448  1.0000     -0.2939 0.3154
Code
#re-analysis of targeting
RItargeting(VQ_progress_data_rasch_recoded_na_omit)

Code
#re-analysis of item hierarchy
itemlabels <- data.frame(itemnr = names(VQ_progress_data_rasch_recoded_na_omit),
                         item = NA)

RIitemHierarchy(VQ_progress_data_rasch_recoded_na_omit)

Code
#re-analysis of threshold locations 
RIitemparams(VQ_progress_data_rasch_recoded_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Item location
VQ.3_recoded -0.78 -0.02 1.00 1.10 NA 0.32
VQ.4_recoded -1.32 -0.86 0.24 NA NA -0.65
VQ.5_recoded -0.19 0.17 0.60 1.42 NA 0.5
VQ.7_recoded -1.35 -0.88 -0.47 0.65 1.47 -0.12
VQ.9_recoded -1.21 -0.50 -0.30 0.47 0.76 -0.16
Note:
Item location is the average of the thresholds for each item.
  • VQ_progress recoded resulted in ordered response categories.
  • item 7 performed slightly better and is kept.
  • item 3 performs worst and bad enough to be removed. It is also the only item with indications of not belonging to the dimension.
  • item 4 also performs bad but is kept until version without item 3 is analyzed.
  • performing Rasch re-analysis with a trimmed recoded VQ_progress without item 3.

5.3 Re-analysis of VQ_progress recoded and trimmed version

Code
#preparing items in VQ_progress recoded and trimmed version
VQ_progress_data_rasch_recoded_trimmed <- bulls_eye_merged_superfinal_data[, c(
  "VQ.4_recoded",
  "VQ.5_recoded",
  "VQ.7_recoded",
  "VQ.9_recoded"
)]

#investigate missing values
RImissing(VQ_progress_data_rasch_recoded_trimmed)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete VQ
VQ_progress_data_rasch_recoded_trimmed_na_omit <- VQ_progress_data_rasch_recoded_trimmed %>%
  select(
    VQ.4_recoded,
    VQ.5_recoded,
    VQ.7_recoded,
    VQ.9_recoded
  ) %>%
  na.omit()

RImissing(VQ_progress_data_rasch_recoded_trimmed_na_omit)
[1] "No missing data."
Code
#tile plot (response category frequency for the items)
RItileplot(VQ_progress_data_rasch_recoded_trimmed_na_omit)

Code
#re-analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI
RItif(na.omit(VQ_progress_data_rasch_recoded_trimmed_na_omit))

Code
RItif(na.omit(VQ_progress_data_rasch_recoded_trimmed_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(VQ_progress_data_rasch_recoded_trimmed_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue sig
overall 17.8 16 0.34
Code
#re-analysis of conditional item fit (assessing unidimensionaliity)
simfit1_VQ_progress_recoded_trimmed <- RIgetfit(VQ_progress_data_rasch_recoded_trimmed_na_omit,
                                       iterations = 100,
                                       cpu = 10) 


RIitemfit(
  VQ_progress_data_rasch_recoded_trimmed_na_omit,
  simfit1_VQ_progress_recoded_trimmed,
  cutoff = c(0.005, 0.995)
)
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
VQ.4_recoded 0.779 [0.834, 1.315] 0.708 [0.771, 1.513] 0.055 0.063 -1.35
VQ.5_recoded 1.152 [0.73, 1.254] 1.125 [0.745, 1.323] no misfit no misfit 0.23
VQ.7_recoded 0.986 [0.788, 1.274] 0.954 [0.812, 1.258] no misfit no misfit -0.66
VQ.9_recoded 1.247 [0.689, 1.31] 1.181 [0.682, 1.375] no misfit no misfit -0.71
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 100 simulated datasets.
Code
#re-analysis of conditional item charachteristic curves
ICCplot(
  as.data.frame(VQ_progress_data_rasch_recoded_trimmed_na_omit),
  itemnumber = 1:4,
  itemdescrip = c(
    "VQ.4_recoded",
    "VQ.5_recoded",
    "VQ.7_recoded",
    "VQ.9_recoded"
  )
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
#re-analysis of PCA of residuals (assessing unidimensionality)
RIpcmPCA(VQ_progress_data_rasch_recoded_trimmed_na_omit)
Eigenvalues Proportion of variance
1.61 43.6%
1.41 34.8%
0.97 21.2%
0.01 0.3%
Code
#re-analysis of item restscores
RIrestscore(VQ_progress_data_rasch_recoded_trimmed_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
VQ.4_recoded 0.73 0.60 0.13 0.059 . -0.77 -1.35
VQ.5_recoded 0.61 0.63 0.02 0.701 0.82 0.23
VQ.7_recoded 0.66 0.63 0.03 0.670 -0.07 -0.66
VQ.9_recoded 0.58 0.64 0.06 0.670 -0.12 -0.71
Code
#re-analysis of loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(VQ_progress_data_rasch_recoded_trimmed_na_omit)

Code
#re-analysis of residual correlations
resid_cor_VQ_progress_recoded_trimmed <- RIgetResidCor(VQ_progress_data_rasch_recoded_trimmed_na_omit,
                                              iterations = 800,
                                              cpu = 10)

RIresidcorr(VQ_progress_data_rasch_recoded_trimmed_na_omit, cutoff = resid_cor_VQ_progress_recoded_trimmed$p99)
VQ.4_recoded VQ.5_recoded VQ.7_recoded VQ.9_recoded
VQ.4_recoded
VQ.5_recoded -0.33
VQ.7_recoded -0.05 -0.19
VQ.9_recoded -0.13 -0.44 -0.52
Note:
Relative cut-off value is 0.038, which is 0.314 above the average correlation (-0.276).
Correlations above the cut-off are highlighted in red text.
Code
#re-analysis of partial gamma coefficients (assessing local independence)
RIpartgamLD(VQ_progress_data_rasch_recoded_trimmed_na_omit) #shows only significant results
Item 1 Item 2 Partial gamma SE Lower CI Upper CI Adjusted p-value (BH)
VQ.4_recoded VQ.7_recoded 0.41 0.143 0.130 0.690 0.049
VQ.4_recoded VQ.9_recoded 0.4 0.133 0.138 0.662 0.033
Code
partgam_LD(as.data.frame(VQ_progress_data_rasch_recoded_trimmed_na_omit)) #shows all output, also non-significant results 
         Item1        Item2   gamma     se pvalue padj.BH sig   lower  upper
1 VQ.4_recoded VQ.5_recoded -0.0405 0.1781 0.8203  1.0000     -0.3895 0.3086
2 VQ.4_recoded VQ.7_recoded  0.4100 0.1428 0.0041  0.0492   *  0.1300 0.6900
3 VQ.4_recoded VQ.9_recoded  0.4000 0.1335 0.0027  0.0328   *  0.1383 0.6617
4 VQ.5_recoded VQ.7_recoded  0.0042 0.1610 0.9791  1.0000     -0.3113 0.3197
5 VQ.5_recoded VQ.9_recoded -0.2245 0.1361 0.0990  1.0000     -0.4912 0.0422
6 VQ.7_recoded VQ.9_recoded -0.0705 0.1426 0.6208  1.0000     -0.3500 0.2089

         Item1        Item2   gamma     se pvalue padj.BH sig   lower   upper
1 VQ.5_recoded VQ.4_recoded -0.3333 0.1689 0.0484  0.5809     -0.6643 -0.0023
2 VQ.7_recoded VQ.4_recoded  0.2577 0.1610 0.1095  1.0000     -0.0579  0.5733
3 VQ.7_recoded VQ.5_recoded  0.1907 0.1500 0.2037  1.0000     -0.1034  0.4847
4 VQ.9_recoded VQ.4_recoded  0.1681 0.1560 0.2812  1.0000     -0.1376  0.4737
5 VQ.9_recoded VQ.5_recoded -0.1349 0.1400 0.3349  1.0000     -0.4093  0.1394
6 VQ.9_recoded VQ.7_recoded -0.2549 0.1427 0.0740  0.8885     -0.5346  0.0248
Code
#re-analysis of targeting
RItargeting(VQ_progress_data_rasch_recoded_trimmed_na_omit)

Code
#re-analysis of item hierarchy
itemlabels <- data.frame(itemnr = names(VQ_progress_data_rasch_recoded_trimmed_na_omit),
                         item = NA)

RIitemHierarchy(VQ_progress_data_rasch_recoded_trimmed_na_omit)

Code
#re-analysis of threshold locations 
RIitemparams(VQ_progress_data_rasch_recoded_trimmed_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Item location
VQ.4_recoded -1.91 -0.98 0.59 NA NA -0.77
VQ.5_recoded -0.41 0.26 1.03 2.39 NA 0.82
VQ.7_recoded -2.08 -1.15 -0.48 0.97 2.38 -0.07
VQ.9_recoded -1.81 -0.73 -0.30 0.73 1.51 -0.12
Note:
Item location is the average of the thresholds for each item.
  • Item 4 performs slightly better but still not optimal. It is however kept in the scale for the validity analysis.
  • The other items works OK but TIF (reliability) is not good at all.

5.4 VQ_progress and BEA validity analysis

Validity between VQ_progress and BEA is investigated by correlational analysis. As both scales has undergone Rasch analysis, participants theta estimates (computed with weighted likelihood estimation - WLE) will be correlated, and since thetas are interval data, Pearson’s correlation is applied.

Code
#complete cases data frame for validity analyses

#removing missing values in relevant items for complete cases analysis
#only BE items and one item from the other scales in the bulls eye study
#is needed (in this case the SubjektivPrestation item) 
#this is because there are no missing values for the other scales (beside some BE items) 
#after the intervention study (Hockeystudien) participants have been removed
#therefore this data frame will work for all validity analyses 

items_validity_analyses <- c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded",
  "SubjektivPrestation"
)

data_complete_cases <- bulls_eye_merged_superfinal_data[complete.cases(bulls_eye_merged_superfinal_data[, items_validity_analyses]), ]

#156 participants to begin with
#34 participants from the intervention study (Hockeystudien) were removed (did not complete any scale beside Bulls Eye)
#8 participants (not from the intervention study) with missing values in Bulls Eye items were removed
#results in a complete cases data frame with 114 participants 

complete_cases_dims <- dim_names(data_complete_cases, "Användarnamn") 
count <- dplyr::count
count(complete_cases_dims)
# A tibble: 1 × 1
      n
  <int>
1   114
Code
#item data frames
BE_validity_analyses <- data_complete_cases[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

VQ_progress_validity_analyses <- data_complete_cases[, c(
  "VQ.4_recoded",
  "VQ.5_recoded",
  "VQ.7_recoded",
  "VQ.9_recoded"
)]

#extracting thetas
BE_thetas <- RIestThetas(BE_validity_analyses)

VQ_progress_thetas <- RIestThetas(VQ_progress_validity_analyses)

#adding thetas to the data frame
data_complete_cases$BE_thetas_WLE <- BE_thetas$WLE 
data_complete_cases$VQ_progress_thetas_WLE <- VQ_progress_thetas$WLE
Code
correlation_BE_VQ_progress <- cor_test(
  data_complete_cases,
  x = "BE_thetas_WLE",
  y = "VQ_progress_thetas_WLE",
  method = "pearson"
)

print(correlation_BE_VQ_progress)
Parameter1    |             Parameter2 |    r |       95% CI | t(112) |         p
---------------------------------------------------------------------------------
BE_thetas_WLE | VQ_progress_thetas_WLE | 0.53 | [0.38, 0.65] |   6.54 | < .001***

Observations: 114
Code
plot(correlation_BE_VQ_progress)

5.4.1 Sensitivity analysis VQ_progress and BEA, time latency

Some participants did not complete BEA and VQ during the same day as intended. Therefore, a sensitivity analysis is conducted.

Code
#creating data frame with only participants that completed Bulls Eye and the other scales during the same day
#since 15 participants did not, 99 participants will remain in the sensitivity analysis
data_complete_cases_less_than_1day <- subset(data_complete_cases, time_latency_survey1_bulls_eye1_days < 1)
Code
#data have previously been prepared for the validity sensitivity analysis
#looking only at participants who completed Bulls Eye and the other scales during the same day (N = 99)
#BE items and BE thetas for these particpants have also been prepared previously

#sensitivity analysis VQ_progress and BE
BE_validity_analyses_less1day <- data_complete_cases_less_than_1day[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

VQ_progress_validity_analyses_less1day <- data_complete_cases_less_than_1day[, c(
  "VQ.4_recoded",
  "VQ.5_recoded",
  "VQ.7_recoded",
  "VQ.9_recoded"
)]

#extracting thetas
BE_thetas_less1day <- RIestThetas(BE_validity_analyses_less1day)

VQ_progress_thetas_less1day <- RIestThetas(VQ_progress_validity_analyses_less1day)

#adding thetas to the data frame
data_complete_cases_less_than_1day$BE_thetas_less1day_WLE <- BE_thetas_less1day$WLE

data_complete_cases_less_than_1day$VQ_progress_thetas_less1day_WLE <- VQ_progress_thetas_less1day$WLE

#analysis only with participants completing Bulls Eye and ADAQ the same day
correlation_BE_VQ_progress_less1day <- cor_test(
  data_complete_cases_less_than_1day,
  x = "BE_thetas_less1day_WLE",
  y = "VQ_progress_thetas_less1day_WLE",
  method = "pearson"
)

print(correlation_BE_VQ_progress_less1day)
Parameter1             |                      Parameter2 |    r |       95% CI
------------------------------------------------------------------------------
BE_thetas_less1day_WLE | VQ_progress_thetas_less1day_WLE | 0.52 | [0.36, 0.65]

Parameter1             | t(97) |         p
------------------------------------------
BE_thetas_less1day_WLE |  6.03 | < .001***

Observations: 99
Code
#we can also check if there is a statistically significant difference 
#by adding same or different day as a predictor 
data_complete_cases_sensitivity_analysis <- data_complete_cases %>% 
  mutate(group = factor(ifelse(time_latency_survey1_bulls_eye1_days < 1,"same_day","different_day")))

m0_VQ_progress <- lm(VQ_progress_thetas_WLE ~ BE_thetas_WLE, data = data_complete_cases_sensitivity_analysis)

m1_VQ_progress <- lm(VQ_progress_thetas_WLE ~ BE_thetas_WLE + group, data = data_complete_cases_sensitivity_analysis)

anova(m0_VQ_progress, m1_VQ_progress)
Analysis of Variance Table

Model 1: VQ_progress_thetas_WLE ~ BE_thetas_WLE
Model 2: VQ_progress_thetas_WLE ~ BE_thetas_WLE + group
  Res.Df    RSS Df Sum of Sq      F  Pr(>F)  
1    112 200.36                              
2    111 194.80  1    5.5636 3.1703 0.07772 .
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Adding whether BE and VQ_progress were completed during the same day or not as a predictor did not significantly alter the results.

Now moving on to analyzing the other part of VQ, VQ_obstruction (to valued living, items 1, 2, 6, 8, 10).

#VQ_obstruction: Rasch analysis and BEA validity investigation

Code
#preparing data frame with VQ_obstruction items
VQ_obstruction_data_rasch <- bulls_eye_merged_superfinal_data[, c("VQ.1", "VQ.2", "VQ.6", "VQ.8", "VQ.10")]

#investigating missing values
RImissing(VQ_obstruction_data_rasch)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete VQ

VQ_obstruction_data_rasch_na_omit <- VQ_obstruction_data_rasch %>%
  select(VQ.1, VQ.2, VQ.6, VQ.8, VQ.10) %>%
  na.omit()

RImissing(VQ_obstruction_data_rasch_na_omit)
[1] "No missing data."
Code
#tile plot (response category frequency for the items)
RItileplot(VQ_obstruction_data_rasch_na_omit)

Code
#analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI 
RItif(na.omit(VQ_obstruction_data_rasch_na_omit))

Code
RItif(na.omit(VQ_obstruction_data_rasch_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(VQ_obstruction_data_rasch_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue sig
overall 32.9 29 0.28
Code
#conditional item fit (assessing unidimensionaliity)
simfit1_VQ_obstruction <- RIgetfit(VQ_obstruction_data_rasch_na_omit,
                                   iterations = 100,
                                   cpu = 10) 

RIitemfit(VQ_obstruction_data_rasch_na_omit, simfit1_VQ_obstruction, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
VQ.1 0.927 [0.806, 1.221] 0.951 [0.79, 1.232] no misfit no misfit -0.28
VQ.2 0.981 [0.743, 1.217] 0.971 [0.765, 1.332] no misfit no misfit -0.42
VQ.6 0.989 [0.773, 1.232] 0.977 [0.725, 1.285] no misfit no misfit -0.53
VQ.8 1.038 [0.812, 1.247] 0.994 [0.779, 1.249] no misfit no misfit -1.10
VQ.10 1.091 [0.81, 1.273] 1.158 [0.798, 1.541] no misfit no misfit -0.51
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 49 simulated datasets.
Code
#conditional item charachteristic curves
ICCplot(
  as.data.frame(VQ_obstruction_data_rasch_na_omit),
  itemnumber = 1:4,
  itemdescrip = c("VQ.1", "VQ.2", "VQ.6", "VQ.8"),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(VQ_obstruction_data_rasch_na_omit),
  itemnumber = 5,
  itemdescrip = c("VQ.10"),
  method = "cut"
)

Code
#PCA of residuals (assessing unidimensionality)
RIpcmPCA(VQ_obstruction_data_rasch_na_omit)
Eigenvalues Proportion of variance
1.69 34.8%
1.27 24.4%
1.07 20.9%
0.97 19.8%
0.01 0.2%
Code
#item restscore
RIrestscore(VQ_obstruction_data_rasch_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
VQ.1 0.47 0.44 0.03 0.991 0.29 -0.28
VQ.2 0.47 0.44 0.03 0.991 0.14 -0.42
VQ.6 0.46 0.44 0.02 0.991 0.04 -0.53
VQ.8 0.40 0.39 0.01 0.991 -0.53 -1.10
VQ.10 0.45 0.45 0.00 0.991 0.06 -0.51
Code
#loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(VQ_obstruction_data_rasch_na_omit)

Code
#residual correlations
resid_cor_VQ_obstruction <- RIgetResidCor(VQ_obstruction_data_rasch_na_omit,
                                          iterations = 800,
                                          cpu = 10)

RIresidcorr(VQ_obstruction_data_rasch_na_omit, cutoff = resid_cor_VQ_obstruction$p99)
VQ.1 VQ.2 VQ.6 VQ.8 VQ.10
VQ.1
VQ.2 -0.01
VQ.6 -0.06 -0.19
VQ.8 -0.28 -0.16 -0.3
VQ.10 -0.35 -0.35 -0.24 0.05
Note:
Relative cut-off value is 0.097, which is 0.286 above the average correlation (-0.189).
Correlations above the cut-off are highlighted in red text.
Code
#partial gamma coefficients (assessing local independence)
RIpartgamLD(VQ_obstruction_data_rasch_na_omit) #shows only significant results
[1] "No statistically significant local dependency found."
Code
partgam_LD(as.data.frame(VQ_obstruction_data_rasch_na_omit)) #shows all output, also non-significant results
   Item1 Item2   gamma     se pvalue padj.BH sig   lower  upper
1   VQ.1  VQ.2  0.1681 0.1483 0.2570       1     -0.1226 0.4587
2   VQ.1  VQ.6  0.1791 0.1224 0.1434       1     -0.0608 0.4190
3   VQ.1  VQ.8 -0.0851 0.1615 0.5982       1     -0.4016 0.2314
4   VQ.1 VQ.10 -0.1184 0.1418 0.4039       1     -0.3963 0.1596
5   VQ.2  VQ.6  0.0424 0.1339 0.7517       1     -0.2201 0.3048
6   VQ.2  VQ.8 -0.0122 0.1570 0.9381       1     -0.3199 0.2955
7   VQ.2 VQ.10  0.0404 0.1374 0.7689       1     -0.2289 0.3096
8   VQ.6  VQ.8 -0.1546 0.1473 0.2938       1     -0.4433 0.1341
9   VQ.6 VQ.10 -0.1811 0.1394 0.1938       1     -0.4543 0.0920
10  VQ.8 VQ.10  0.2162 0.1543 0.1612       1     -0.0863 0.5187

   Item1 Item2   gamma     se pvalue padj.BH sig   lower   upper
1   VQ.2  VQ.1  0.2544 0.1305 0.0512  1.0000     -0.0014  0.5101
2   VQ.6  VQ.1  0.2297 0.1243 0.0646  1.0000     -0.0139  0.4734
3   VQ.6  VQ.2 -0.0123 0.1336 0.9264  1.0000     -0.2743  0.2496
4   VQ.8  VQ.1 -0.3950 0.1319 0.0028  0.0552   . -0.6536 -0.1364
5   VQ.8  VQ.2 -0.0974 0.1480 0.5102  1.0000     -0.3875  0.1926
6   VQ.8  VQ.6 -0.1142 0.1493 0.4444  1.0000     -0.4067  0.1784
7  VQ.10  VQ.1 -0.2483 0.1308 0.0578  1.0000     -0.5047  0.0082
8  VQ.10  VQ.2 -0.0386 0.1484 0.7946  1.0000     -0.3294  0.2521
9  VQ.10  VQ.6 -0.1613 0.1430 0.2595  1.0000     -0.4416  0.1191
10 VQ.10  VQ.8  0.1919 0.1441 0.1828  1.0000     -0.0904  0.4743
Code
#targeting
RItargeting(VQ_obstruction_data_rasch_na_omit)

Code
#item hierarchy
itemlabels <- data.frame(itemnr = names(VQ_obstruction_data_rasch_na_omit),
                         item = NA)

RIitemHierarchy(VQ_obstruction_data_rasch_na_omit)

Code
#threshold locations 
RIitemparams(VQ_obstruction_data_rasch_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Threshold 5 Threshold 6 Item location
VQ.1 -0.06 -0.75 -0.16 0.74 0.86 1.10 0.29
VQ.2 0.41 -1.85 0.05 0.58 1.24 0.42 0.14
VQ.6 -1.30 -0.60 0.39 0.54 0.20 1.01 0.04
VQ.8 -0.38 -2.13 0.02 -0.87 0.13 0.05 -0.53
VQ.10 -1.10 -0.28 0.06 0.70 0.24 0.74 0.06
Note:
Item location is the average of the thresholds for each item.
Code
#analysis of response categories
RIitemCats(
  VQ_obstruction_data_rasch_na_omit,
  items = "all",
  xlims = c(-5, 5),
  legend = "left"
)

5.5 Recoding response categories for VQ_obstruction

Based on the analysis of the response categories and threshold locations, scale steps were merged and recoded so that each kept scale step had the highest response probability at some point on the person location scale, and also so that all threshold locations were in the right order for all items. Recoded items were given new names with the abbrevation “_recoded”.

Code
#add recoded VQ_obstruction items

bulls_eye_merged_superfinal_data <- bulls_eye_merged_superfinal_data %>%
  mutate(
    # Recoding VQ.1 (0+1+2)
    VQ.1_recoded = recode(VQ.1,
                            `0` = 0, `1` = 0, `2` = 0, 
                            `3` = 1,
                            `4` = 2,
                            `5` = 3,
                            `6` = 4),
    # Recoding VQ.2 (0+1; 4+5)
    VQ.2_recoded = recode(VQ.2,
                            `0` = 0, `1` = 0, 
                            `2` = 1,
                            `3` = 2, 
                            `4` = 3, `5` = 3,
                            `6` = 4),
    # Recoding VQ.6 (0+1; 3+4)
    VQ.6_recoded = recode(VQ.6,
                           `0` = 0, `1` = 0,
                           `2` = 1, 
                           `3` = 2, `4` = 2,
                           `5` = 3,
                           `6` = 4),
    # Recoding VQ.8 (0+1+2+3; 4+5)
    VQ.8_recoded = recode(VQ.8,
                           `0` = 0, `1` = 0, `2` = 0, `3` = 0, 
                           `4` = 1, `5` = 1,
                           `6` = 2),
    # Recoding VQ.10 (0+1, 2+3; 4+5)
    VQ.10_recoded = recode(VQ.10,
                           `0` = 0, `1` = 0,
                           `2` = 1, `3` = 1,
                           `4` = 2, `5` = 2,
                           `6` = 3))

5.6 Re-analysis of VQ_obstruction recoded version

Code
VQ_obstruction_data_rasch_recoded <- bulls_eye_merged_superfinal_data[, c(
  "VQ.1_recoded",
  "VQ.2_recoded",
  "VQ.6_recoded",
  "VQ.8_recoded",
  "VQ.10_recoded"
)]

#investigating missing values
RImissing(VQ_obstruction_data_rasch_recoded)

Code
#removing missing values 
#in practice this means removing participants from the intervention study (Hockeystudien) that did not complete VQ
VQ_obstruction_data_rasch_recoded_na_omit <- VQ_obstruction_data_rasch_recoded %>%
  select(
    VQ.1_recoded,
    VQ.2_recoded,
    VQ.6_recoded,
    VQ.8_recoded,
    VQ.10_recoded
  ) %>%
  na.omit()

RImissing(VQ_obstruction_data_rasch_recoded_na_omit)
[1] "No missing data."
Code
#re-analysis of response categories using VQ_obstruction recoded
RIitemCats(
  VQ_obstruction_data_rasch_recoded_na_omit,
  items = "all",
  xlims = c(-5, 5),
  legend = "left"
)

Code
#tile plot (response category frequency for the items)
RItileplot(VQ_obstruction_data_rasch_recoded_na_omit)

Code
#re-analysis of reliability with Test Information Function (TIF)
#first output shows without sample PSI and second output with sample PSI 
RItif(na.omit(VQ_obstruction_data_rasch_recoded_na_omit))

Code
RItif(na.omit(VQ_obstruction_data_rasch_recoded_na_omit), samplePSI = TRUE)

Code
#conditional likelihood ratio tests (CLR)
#a test of homogeneity that assesses whether item parameters are homogeneous across all levels of theta 
clr_tests(VQ_obstruction_data_rasch_recoded_na_omit, model = "PCM")

Conditional Likelihood Ratio Tests:
        clr  df pvalue sig
overall 21.9 16 0.15
Code
#re-analysis of conditional item fit (assessing unidimensionaliity)
simfit1_VQ_obstruction_recoded <- RIgetfit(VQ_obstruction_data_rasch_recoded_na_omit,
                                           iterations = 100,
                                           cpu = 10)

RIitemfit(VQ_obstruction_data_rasch_recoded_na_omit, simfit1_VQ_obstruction_recoded, cutoff = c(0.005, 0.995))
Item InfitMSQ Infit thresholds OutfitMSQ Outfit thresholds Infit diff Outfit diff Relative location
VQ.1_recoded 1.091 [0.816, 1.231] 1.01 [0.806, 1.32] no misfit no misfit 0.23
VQ.2_recoded 0.997 [0.746, 1.309] 0.995 [0.745, 1.399] no misfit no misfit -0.47
VQ.6_recoded 0.968 [0.756, 1.258] 0.945 [0.755, 1.223] no misfit no misfit -0.27
VQ.8_recoded 1.039 [0.765, 1.237] 1.003 [0.788, 1.271] no misfit no misfit -0.79
VQ.10_recoded 1.016 [0.769, 1.169] 1.05 [0.747, 1.23] no misfit no misfit -0.31
Note:
MSQ values based on conditional calculations (n = 122 complete cases).
Simulation based thresholds from 100 simulated datasets.
Code
#re-analysis of conditional item charachteristic curves
ICCplot(
  as.data.frame(VQ_obstruction_data_rasch_recoded_na_omit),
  itemnumber = 1:4,
  itemdescrip = c(
    "VQ.1_recoded",
    "VQ.2_recoded",
    "VQ.6_recoded",
    "VQ.8_recoded"
  ),
  method = "cut"
)

[1] "Please press Zoom on the Plots window to see the plot"
Code
ICCplot(
  as.data.frame(VQ_obstruction_data_rasch_recoded_na_omit),
  itemnumber = 5,
  itemdescrip = c("VQ.10_recoded"),
  method = "cut"
)

Code
#re-analyzing PCA of residuals (assessing unidimensionality)
RIpcmPCA(VQ_obstruction_data_rasch_recoded_na_omit)
Eigenvalues Proportion of variance
1.61 33.4%
1.32 25.1%
1.17 22.4%
0.88 19%
0.01 0.1%
Code
#re-analysis of item restscores
RIrestscore(VQ_obstruction_data_rasch_recoded_na_omit)
Item Observed value Model expected value Absolute difference Adjusted p-value (BH) Statistical significance level Location Relative location
VQ.1_recoded 0.50 0.49 0.01 0.926 0.50 0.23
VQ.2_recoded 0.50 0.46 0.04 0.926 -0.20 -0.47
VQ.6_recoded 0.48 0.47 0.01 0.926 0.00 -0.27
VQ.8_recoded 0.37 0.40 0.03 0.926 -0.52 -0.79
VQ.10_recoded 0.46 0.44 0.02 0.926 -0.05 -0.31
Code
#re-analysis of loadings on first residual contrast (assessing clusters in data or multidimensionality)
RIloadLoc(VQ_obstruction_data_rasch_recoded_na_omit)

Code
#re-analysis of residual correlations
resid_cor_VQ_obstruction_recoded <- RIgetResidCor(VQ_obstruction_data_rasch_recoded_na_omit,
                                                  iterations = 800,
                                                  cpu = 10)

RIresidcorr(VQ_obstruction_data_rasch_recoded_na_omit, cutoff = resid_cor_VQ_obstruction_recoded$p99)
VQ.1_recoded VQ.2_recoded VQ.6_recoded VQ.8_recoded VQ.10_recoded
VQ.1_recoded
VQ.2_recoded -0.15
VQ.6_recoded -0.16 -0.25
VQ.8_recoded -0.35 -0.14 -0.15
VQ.10_recoded -0.29 -0.22 -0.23 0.13
Note:
Relative cut-off value is 0.147, which is 0.33 above the average correlation (-0.182).
Correlations above the cut-off are highlighted in red text.
Code
#re-analysis of partial gamma coefficients (assessing local independence)
RIpartgamLD(VQ_obstruction_data_rasch_recoded_na_omit) #shows only significant results
[1] "No statistically significant local dependency found."
Code
partgam_LD(as.data.frame(VQ_obstruction_data_rasch_recoded_na_omit)) #shows all output, also non-significant results
          Item1         Item2   gamma     se pvalue padj.BH sig   lower  upper
1  VQ.1_recoded  VQ.2_recoded  0.1521 0.1423 0.2852       1     -0.1268 0.4311
2  VQ.1_recoded  VQ.6_recoded  0.1754 0.1309 0.1803       1     -0.0811 0.4319
3  VQ.1_recoded  VQ.8_recoded -0.2906 0.1659 0.0797       1     -0.6157 0.0345
4  VQ.1_recoded VQ.10_recoded -0.1321 0.1571 0.4005       1     -0.4400 0.1758
5  VQ.2_recoded  VQ.6_recoded -0.0168 0.1382 0.9031       1     -0.2878 0.2541
6  VQ.2_recoded  VQ.8_recoded -0.0806 0.1669 0.6293       1     -0.4077 0.2466
7  VQ.2_recoded VQ.10_recoded  0.0000 0.1464 1.0000       1     -0.2869 0.2869
8  VQ.6_recoded  VQ.8_recoded -0.0132 0.1632 0.9355       1     -0.3331 0.3066
9  VQ.6_recoded VQ.10_recoded -0.1055 0.1549 0.4960       1     -0.4090 0.1981
10 VQ.8_recoded VQ.10_recoded  0.2750 0.1681 0.1019       1     -0.0545 0.6045

           Item1        Item2   gamma     se pvalue padj.BH sig   lower   upper
1   VQ.2_recoded VQ.1_recoded  0.1180 0.1310 0.3679  1.0000     -0.1388  0.3748
2   VQ.6_recoded VQ.1_recoded  0.1599 0.1241 0.1977  1.0000     -0.0834  0.4032
3   VQ.6_recoded VQ.2_recoded -0.0789 0.1373 0.5654  1.0000     -0.3481  0.1903
4   VQ.8_recoded VQ.1_recoded -0.3433 0.1473 0.0198  0.3953     -0.6320 -0.0546
5   VQ.8_recoded VQ.2_recoded -0.0168 0.1709 0.9217  1.0000     -0.3518  0.3182
6   VQ.8_recoded VQ.6_recoded -0.0866 0.1577 0.5828  1.0000     -0.3957  0.2224
7  VQ.10_recoded VQ.1_recoded -0.1026 0.1452 0.4799  1.0000     -0.3871  0.1820
8  VQ.10_recoded VQ.2_recoded -0.0609 0.1648 0.7116  1.0000     -0.3839  0.2621
9  VQ.10_recoded VQ.6_recoded -0.1447 0.1508 0.3372  1.0000     -0.4403  0.1509
10 VQ.10_recoded VQ.8_recoded  0.3488 0.1661 0.0357  0.7136      0.0233  0.6743
Code
#re-analysis of targeting
RItargeting(VQ_obstruction_data_rasch_recoded_na_omit)

Code
#re-analysis of item hierarchy
itemlabels <- data.frame(itemnr = names(VQ_obstruction_data_rasch_recoded_na_omit),
                         item = NA)

RIitemHierarchy(VQ_obstruction_data_rasch_recoded_na_omit)

Code
#threshold locations 
RIitemparams(VQ_obstruction_data_rasch_recoded_na_omit)
Threshold 1 Threshold 2 Threshold 3 Threshold 4 Item location
VQ.1_recoded -0.08 0.38 0.68 1.02 0.5
VQ.2_recoded -1.47 -0.36 -0.11 1.13 -0.2
VQ.6_recoded -0.83 -0.60 0.58 0.83 0
VQ.8_recoded -1.29 0.25 NA NA -0.52
VQ.10_recoded -1.14 0.06 0.94 NA -0.05
Note:
Item location is the average of the thresholds for each item.

Response categories in VQ_obstruction is now ordered. All items in VQ_obstruction works fine enough for all to be kept. Moving on to the validity analysis.

5.7 VQ_obstruction and BEA validity analysis

Validity between VQ_obstruction and BEA is investigated by correlational analysis. As both scales has undergone Rasch analysis, participants theta estimates (computed with weighted likelihood estimation - WLE) will be correlated, and since thetas are interval data, Pearson’s correlation is applied.

Code
#complete cases data frame for validity analyses

#removing missing values in relevant items for complete cases analysis
#only BE items and one item from the other scales in the bulls eye study
#is needed (in this case the SubjektivPrestation item) 
#this is because there are no missing values for the other scales (beside some BE items) 
#after the intervention study (Hockeystudien) participants have been removed
#therefore this data frame will work for all validity analyses 

items_validity_analyses <- c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded",
  "SubjektivPrestation"
)

data_complete_cases <- bulls_eye_merged_superfinal_data[complete.cases(bulls_eye_merged_superfinal_data[, items_validity_analyses]), ]

#156 participants to begin with
#34 participants from the intervention study (Hockeystudien) were removed (did not complete any scale beside Bulls Eye)
#8 participants (not from the intervention study) with missing values in Bulls Eye items were removed
#results in a complete cases data frame with 114 participants 

complete_cases_dims <- dim_names(data_complete_cases, "Användarnamn") 
count <- dplyr::count
count(complete_cases_dims)
# A tibble: 1 × 1
      n
  <int>
1   114
Code
#item data frames
BE_validity_analyses <- data_complete_cases[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

VQ_obstruction_validity_analyses <- data_complete_cases[, c(
  "VQ.1_recoded",
  "VQ.2_recoded",
  "VQ.6_recoded",
  "VQ.8_recoded",
  "VQ.10_recoded"
)]

#extracting thetas
BE_thetas <- RIestThetas(BE_validity_analyses)

VQ_obstruction_thetas <- RIestThetas(VQ_obstruction_validity_analyses)

#adding thetas to the data frame
data_complete_cases$BE_thetas_WLE <- BE_thetas$WLE 
data_complete_cases$VQ_obstruction_thetas_WLE <- VQ_obstruction_thetas$WLE
Code
correlation_BE_VQ_obstruction <- cor_test(
  data_complete_cases,
  x = "BE_thetas_WLE",
  y = "VQ_obstruction_thetas_WLE",
  method = "pearson"
)

print(correlation_BE_VQ_obstruction)
Parameter1    |                Parameter2 |    r |       95% CI | t(112) |         p
------------------------------------------------------------------------------------
BE_thetas_WLE | VQ_obstruction_thetas_WLE | 0.48 | [0.32, 0.61] |   5.71 | < .001***

Observations: 114
Code
plot(correlation_BE_VQ_obstruction)

5.7.1 Sensitivity analysis VQ_obstruction and BEA, time latency

Some participants did not complete BEA and VQ during the same day as intended. Therefore, a sensitivity analysis is conducted.

Code
#creating data frame with only participants that completed Bulls Eye and the other scales during the same day
#since 15 participants did not, 99 participants will remain in the sensitivity analysis
data_complete_cases_less_than_1day <- subset(data_complete_cases, time_latency_survey1_bulls_eye1_days < 1)
Code
#data have previously been prepared for the validity sensitivity analysis
#looking only at participants who completed Bulls Eye and the other scales during the same day (N = 99)
#BE items and BE thetas for these particpants have also been prepared previously

#sensitivity analysis VQ_obstruction and BE
BE_validity_analyses_less1day <- data_complete_cases_less_than_1day[, c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded"
)]

VQ_obstruction_validity_analyses_less1day <- data_complete_cases_less_than_1day[, c(
  "VQ.1_recoded",
  "VQ.2_recoded",
  "VQ.6_recoded",
  "VQ.8_recoded",
  "VQ.10_recoded"
)]

#extracting thetas
BE_thetas_less1day <- RIestThetas(BE_validity_analyses_less1day)

VQ_obstruction_thetas_less1day <- RIestThetas(VQ_obstruction_validity_analyses_less1day)

#adding thetas to the data frame
data_complete_cases_less_than_1day$BE_thetas_less1day_WLE <- BE_thetas_less1day$WLE

data_complete_cases_less_than_1day$VQ_obstruction_thetas_less1day_WLE <- VQ_obstruction_thetas_less1day$WLE

#analysis only with participants completing Bulls Eye and ADAQ the same day
correlation_BE_VQ_obstruction_less1day <- cor_test(
  data_complete_cases_less_than_1day,
  x = "BE_thetas_less1day_WLE",
  y = "VQ_obstruction_thetas_less1day_WLE",
  method = "pearson"
)

print(correlation_BE_VQ_obstruction_less1day)
Parameter1             |                         Parameter2 |    r
------------------------------------------------------------------
BE_thetas_less1day_WLE | VQ_obstruction_thetas_less1day_WLE | 0.45

Parameter1             |       95% CI | t(97) |         p
---------------------------------------------------------
BE_thetas_less1day_WLE | [0.27, 0.59] |  4.91 | < .001***

Observations: 99
Code
#we can also check if there is a statistically significant difference 
#by adding same or different day as a predictor 
data_complete_cases_sensitivity_analysis <- data_complete_cases %>% 
  mutate(group = factor(ifelse(time_latency_survey1_bulls_eye1_days < 1,"same_day","different_day")))

m0_VQ_obstruction <- lm(VQ_obstruction_thetas_WLE ~ BE_thetas_WLE, data = data_complete_cases_sensitivity_analysis)

m1_VQ_obstruction <- lm(VQ_obstruction_thetas_WLE ~ BE_thetas_WLE + group, data = data_complete_cases_sensitivity_analysis)

anova(m0_VQ_obstruction, m1_VQ_obstruction)
Analysis of Variance Table

Model 1: VQ_obstruction_thetas_WLE ~ BE_thetas_WLE
Model 2: VQ_obstruction_thetas_WLE ~ BE_thetas_WLE + group
  Res.Df    RSS Df Sum of Sq      F Pr(>F)
1    112 83.076                           
2    111 83.042  1  0.034302 0.0459 0.8308

Adding whether BE and VQ_obstruction were completed during the same day or not as a predictor did not significantly alter the results.

6 Sport Motivation Scale-2 and BEA validity analysis

Validity between Sport Motivation Scale-2 (SMS-2) and BEA is investigated by correlational analysis. SMS-2 has not undergone Rasch analysis, thus participants sum scores will be correlated. Since sum scores are ordinal data, Spearman’s rank correlation is applied.

Code
#complete cases data frame for validity analyses

#removing missing values in relevant items for complete cases analysis
#only BE items and one item from the other scales in the bulls eye study
#is needed (in this case the SubjektivPrestation item) 
#this is because there are no missing values for the other scales (beside some BE items) 
#after the intervention study (Hockeystudien) participants have been removed
#therefore this data frame will work for all validity analyses 

items_validity_analyses <- c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded",
  "SubjektivPrestation"
)

data_complete_cases <- bulls_eye_merged_superfinal_data[complete.cases(bulls_eye_merged_superfinal_data[, items_validity_analyses]), ]

#156 participants to begin with
#34 participants from the intervention study (Hockeystudien) were removed (did not complete any scale beside Bulls Eye)
#8 participants (not from the intervention study) with missing values in Bulls Eye items were removed
#results in a complete cases data frame with 114 participants 

complete_cases_dims <- dim_names(data_complete_cases, "Användarnamn") 
count <- dplyr::count
count(complete_cases_dims)
# A tibble: 1 × 1
      n
  <int>
1   114
Code
Code
#item preparation and plot of sum scores 
SMS.2_items <- data_complete_cases[, c(
  "SMS-2.1",
  "SMS-2.2",
  "SMS-2.3",
  "SMS-2.4",
  "SMS-2.5",
  "SMS-2.6",
  "SMS-2.7",
  "SMS-2.8",
  "SMS-2.9",
  "SMS-2.10",
  "SMS-2.11",
  "SMS-2.12",
  "SMS-2.13",
  "SMS-2.14",
  "SMS-2.15",
  "SMS-2.16",
  "SMS-2.17",
  "SMS-2.18"
)]


BE_items <- data_complete_cases[, c(
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life"
)]

BE_items_names <- names(BE_items)
SMS.2_items_names <- names(SMS.2_items)

#computing sum scores for each scale
data_complete_cases_SMS.2_BE <- data_complete_cases %>%
  mutate(SMS_2_sum = rowSums(select(., all_of(SMS.2_items_names)), na.rm = TRUE),
         BE_sum = rowSums(select(., all_of(BE_items_names)), na.rm = TRUE))

#plot of sum scores (raw data, not ranked)
ggplot(data_complete_cases_SMS.2_BE, aes(x = SMS_2_sum, y = BE_sum)) +
  geom_point() +
  labs(
    title = "Scatterplot of SMS-2 vs BE Sum Scores",
    x = "SMS-2 Sum Score",
    y = "BE Sum Score"
  ) +
  theme_minimal()

Code
#compute Spearman's rank correlation
correlation_SMS.2_BE <- correlation(data_complete_cases_SMS.2_BE[, c("SMS_2_sum", "BE_sum")], method = "spearman")
print(correlation_SMS.2_BE)
# Correlation Matrix (spearman-method)

Parameter1 | Parameter2 |  rho |        95% CI |        S |     p
-----------------------------------------------------------------
SMS_2_sum  |     BE_sum | 0.11 | [-0.08, 0.30] | 2.19e+05 | 0.231

p-value adjustment method: Holm (1979)
Observations: 114

6.1 Sensitivity analysis SMS-2 and BEA, time latency

Some participants did not complete BEA and SMS-2 during the same day as intended. Therefore, a sensitivity analysis is conducted.

Code
#creating data frame with only participants that completed Bulls Eye and the other scales during the same day
#since 15 participants did not, 99 participants will remain in the sensitivity analysis
data_complete_cases_less_than_1day <- subset(data_complete_cases, time_latency_survey1_bulls_eye1_days < 1)
Code
#sensitivity analysis correlation SMS-2 and BEA
SMS.2_items_less1day <- data_complete_cases_less_than_1day[, c(
  "SMS-2.1",
  "SMS-2.2",
  "SMS-2.3",
  "SMS-2.4",
  "SMS-2.5",
  "SMS-2.6",
  "SMS-2.7",
  "SMS-2.8",
  "SMS-2.9",
  "SMS-2.10",
  "SMS-2.11",
  "SMS-2.12",
  "SMS-2.13",
  "SMS-2.14",
  "SMS-2.15",
  "SMS-2.16",
  "SMS-2.17",
  "SMS-2.18"
)]


BE_items_less1day <- data_complete_cases_less_than_1day[, c(
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life"
)]

BE_items_less1day_names <- names(BE_items_less1day)
SMS.2_items_less1day_names <- names(SMS.2_items_less1day)


#computing sum scores for each scale
data_complete_cases_SMS.2_BE_less1day <- data_complete_cases_less_than_1day %>%
  mutate(SMS_2_sum = rowSums(select(., all_of(
    SMS.2_items_less1day_names
  )), na.rm = TRUE),
  BE_sum = rowSums(select(., all_of(
    BE_items_less1day_names
  )), na.rm = TRUE))

#compute Spearman's rank correlation
correlation_SMS.2_BE_less1day <- correlation(data_complete_cases_SMS.2_BE_less1day[, c("SMS_2_sum", "BE_sum")], method = "spearman")

print(correlation_SMS.2_BE_less1day)
# Correlation Matrix (spearman-method)

Parameter1 | Parameter2 |  rho |        95% CI |        S |     p
-----------------------------------------------------------------
SMS_2_sum  |     BE_sum | 0.11 | [-0.10, 0.30] | 1.45e+05 | 0.297

p-value adjustment method: Holm (1979)
Observations: 99

Result from sensitivity analysis correlation is very close to the original correlation.

7 Competitive State Anxiety Inverntory-2 Revised and BEA validity analysis

Validity between Competitive State Anxiety Inverntory-2 Revised (CSAI-2R) and BEA is investigated by correlational analysis. CSAI-2R has not undergone Rasch analysis, thus participants sum scores will be correlated. Since sum scores are ordinal data, Spearman’s rank correlation is applied.

Code
#complete cases data frame for validity analyses

#removing missing values in relevant items for complete cases analysis
#only BE items and one item from the other scales in the bulls eye study
#is needed (in this case the SubjektivPrestation item) 
#this is because there are no missing values for the other scales (beside some BE items) 
#after the intervention study (Hockeystudien) participants have been removed
#therefore this data frame will work for all validity analyses 

items_validity_analyses <- c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded",
  "SubjektivPrestation"
)

data_complete_cases <- bulls_eye_merged_superfinal_data[complete.cases(bulls_eye_merged_superfinal_data[, items_validity_analyses]), ]

#156 participants to begin with
#34 participants from the intervention study (Hockeystudien) were removed (did not complete any scale beside Bulls Eye)
#8 participants (not from the intervention study) with missing values in Bulls Eye items were removed
#results in a complete cases data frame with 114 participants 

complete_cases_dims <- dim_names(data_complete_cases, "Användarnamn") 
count <- dplyr::count
count(complete_cases_dims)
# A tibble: 1 × 1
      n
  <int>
1   114
Code
#item preparation and plot of sum scores
CSAI2R_items <- data_complete_cases[, c(
  "CSAI-2R.1",
  "CSAI-2R.2",
  "CSAI-2R.3",
  "CSAI-2R.4",
  "CSAI-2R.5",
  "CSAI-2R.6",
  "CSAI-2R.7",
  "CSAI-2R.8",
  "CSAI-2R.9",
  "CSAI-2R.10",
  "CSAI-2R.11",
  "CSAI-2R.12",
  "CSAI-2R.13",
  "CSAI-2R.14",
  "CSAI-2R.15",
  "CSAI-2R.16",
  "CSAI-2R.17"
)]

BE_items <- data_complete_cases[, c(
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life"
)]

BE_items_names <- names(BE_items)
CSAI2R_items_names <- names(CSAI2R_items)

#computing sum scores for each scale
data_complete_cases_CSAI2R_BE <- data_complete_cases %>%
  mutate(CSAI2R_sum = rowSums(select(., all_of(CSAI2R_items_names)), na.rm = TRUE),
         BE_sum = rowSums(select(., all_of(BE_items_names)), na.rm = TRUE))

#plot of sum scores (raw data, not ranked)
ggplot(data_complete_cases_CSAI2R_BE, aes(x = CSAI2R_sum, y = BE_sum)) +
  geom_point() +
  labs(
    title = "Scatterplot of Performance rating vs BE Sum Scores",
    x = "CSAI-2R Sum Score",
    y = "BE Sum Score"
  ) +
  theme_minimal()

Code
#compute Spearman's rank correlation
correlation_CSAI2R_BE <- correlation(data_complete_cases_CSAI2R_BE[, c("CSAI2R_sum", "BE_sum")], method = "spearman")

print(correlation_CSAI2R_BE)
# Correlation Matrix (spearman-method)

Parameter1 | Parameter2 |   rho |         95% CI |        S |         p
-----------------------------------------------------------------------
CSAI2R_sum |     BE_sum | -0.41 | [-0.55, -0.23] | 3.47e+05 | < .001***

p-value adjustment method: Holm (1979)
Observations: 114

7.1 Sensitivity analysis CSAI-2R and BEA, time latency

Some participants did not complete BEA and CSAI-2R during the same day as intended. Therefore, a sensitivity analysis is conducted.

Code
#creating data frame with only participants that completed Bulls Eye and the other scales during the same day
#since 15 participants did not, 99 participants will remain in the sensitivity analysis
data_complete_cases_less_than_1day <- subset(data_complete_cases, time_latency_survey1_bulls_eye1_days < 1)
Code
#sensitivity analysis correlation BEA and CSAI-2R
CSAI2R_items_less1day <- data_complete_cases_less_than_1day[, c(
  "CSAI-2R.1",
  "CSAI-2R.2",
  "CSAI-2R.3",
  "CSAI-2R.4",
  "CSAI-2R.5",
  "CSAI-2R.6",
  "CSAI-2R.7",
  "CSAI-2R.8",
  "CSAI-2R.9",
  "CSAI-2R.10",
  "CSAI-2R.11",
  "CSAI-2R.12",
  "CSAI-2R.13",
  "CSAI-2R.14",
  "CSAI-2R.15",
  "CSAI-2R.16",
  "CSAI-2R.17"
)]

BE_items_less1day <- data_complete_cases_less_than_1day[, c(
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life"
)]

BE_items_less1day_names <- names(BE_items_less1day)
CSAI2R_items_less1day_names <- names(CSAI2R_items_less1day)

#computing sum scores for each scale
data_complete_cases_CSAI2R_BE_less1day <- data_complete_cases_less_than_1day %>%
  mutate(
    CSAI2R_sum = rowSums(select(., all_of(
      CSAI2R_items_less1day_names
    )), na.rm = TRUE),
    BE_sum = rowSums(select(., all_of(
      BE_items_less1day_names
    )), na.rm = TRUE)
  )

#compute Spearman's rank correlation
correlation_CSAI2R_BE_less1day <- correlation(data_complete_cases_CSAI2R_BE_less1day[, c("CSAI2R_sum", "BE_sum")], method = "spearman")

print(correlation_CSAI2R_BE_less1day)
# Correlation Matrix (spearman-method)

Parameter1 | Parameter2 |   rho |         95% CI |        S |         p
-----------------------------------------------------------------------
CSAI2R_sum |     BE_sum | -0.48 | [-0.62, -0.31] | 2.40e+05 | < .001***

p-value adjustment method: Holm (1979)
Observations: 99

Result from sensitivity analysis correlation is somewhat close to the original correlation, at least interpretation wise.

8 Subjective performance rating and BEA validity analysis

Validity between the subjective performance rating and BEA is investigated by correlational analysis. The subjective performance rating has not undergone Rasch analysis (since it is only one item), thus participants sum scores will be correlated. Since sum scores are ordinal data, Spearman’s rank correlation is applied.

Code
#complete cases data frame for validity analyses

#removing missing values in relevant items for complete cases analysis
#only BE items and one item from the other scales in the bulls eye study
#is needed (in this case the SubjektivPrestation item) 
#this is because there are no missing values for the other scales (beside some BE items) 
#after the intervention study (Hockeystudien) participants have been removed
#therefore this data frame will work for all validity analyses 

items_validity_analyses <- c(
  "BE_Competition_recoded",
  "BE_Training_recoded",
  "BE_PreparationRecovery_recoded",
  "BE_Life_recoded",
  "BE_Obstacles_Sport_recoded",
  "BE_Obstacles_Life_recoded",
  "SubjektivPrestation"
)

data_complete_cases <- bulls_eye_merged_superfinal_data[complete.cases(bulls_eye_merged_superfinal_data[, items_validity_analyses]), ]

#156 participants to begin with
#34 participants from the intervention study (Hockeystudien) were removed (did not complete any scale beside Bulls Eye)
#8 participants (not from the intervention study) with missing values in Bulls Eye items were removed
#results in a complete cases data frame with 114 participants 

complete_cases_dims <- dim_names(data_complete_cases, "Användarnamn") 
count <- dplyr::count
count(complete_cases_dims)
# A tibble: 1 × 1
      n
  <int>
1   114
Code
#item preparation and plot of sum scores
performance_rating <- data_complete_cases[, c("SubjektivPrestation")]

BE_items <- data_complete_cases[, c(
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life"
)]

BE_items_names <- names(BE_items)
performance_rating_names <- names(performance_rating)

#computing sum scores 
data_complete_cases_performance_rating_BE <- data_complete_cases %>%
  mutate(
    performance_rating_sum = rowSums(select(., all_of(
      performance_rating_names
    )), na.rm = TRUE),
    BE_sum = rowSums(select(., all_of(BE_items_names)), na.rm = TRUE)
  )

#plot of sum scores (raw data, not ranked)
ggplot(data_complete_cases_performance_rating_BE, aes(x = performance_rating_sum, y = BE_sum)) +
  geom_point() +
  labs(
    title = "Scatterplot of Performance rating vs BE Sum Scores",
    x = "Performance rating",
    y = "BE Sum Score"
  ) +
  theme_minimal()

Code
#compute Spearman's rank correlation
correlation_performance_rating_BE <- correlation(data_complete_cases_performance_rating_BE[, c("performance_rating_sum", "BE_sum")], method = "spearman")

print(correlation_performance_rating_BE)
# Correlation Matrix (spearman-method)

Parameter1             | Parameter2 |  rho |       95% CI |        S |       p
------------------------------------------------------------------------------
performance_rating_sum |     BE_sum | 0.27 | [0.09, 0.44] | 1.79e+05 | 0.003**

p-value adjustment method: Holm (1979)
Observations: 114

8.1 Sensitivity analysis Subjective performance rating and BEA, time latency

Some participants did not complete BEA and the subjective performance rating during the same day as intended. Therefore, a sensitivity analysis is conducted.

Code
#creating data frame with only participants that completed Bulls Eye and the other scales during the same day
#since 15 participants did not, 99 participants will remain in the sensitivity analysis
data_complete_cases_less_than_1day <- subset(data_complete_cases, time_latency_survey1_bulls_eye1_days < 1)
Code
#sensitivity correlation analysis BEA and subjective performance rating
performance_rating_less1day <- data_complete_cases_less_than_1day[, c("SubjektivPrestation")]

BE_items_less1day <- data_complete_cases_less_than_1day[, c(
  "BE_Competition",
  "BE_Training",
  "BE_PreparationRecovery",
  "BE_Life",
  "BE_Obstacles_Sport",
  "BE_Obstacles_Life"
)]

BE_items_less1day_names <- names(BE_items_less1day)
performance_rating_less1day_names <- names(performance_rating_less1day)

#computing sum scores
data_complete_cases_performance_rating_BE_less1day <- data_complete_cases_less_than_1day %>%
  mutate(
    performance_rating_sum = rowSums(select(
      ., all_of(performance_rating_less1day_names)
    ), na.rm = TRUE),
    BE_sum = rowSums(select(., all_of(
      BE_items_less1day_names
    )), na.rm = TRUE)
  )

#compute Spearman's rank correlation
correlation_performance_rating_BE_less1day <- correlation(data_complete_cases_performance_rating_BE_less1day[, c("performance_rating_sum", "BE_sum")], method = "spearman")

print(correlation_performance_rating_BE_less1day)
# Correlation Matrix (spearman-method)

Parameter1             | Parameter2 |  rho |       95% CI |        S |       p
------------------------------------------------------------------------------
performance_rating_sum |     BE_sum | 0.29 | [0.09, 0.46] | 1.15e+05 | 0.004**

p-value adjustment method: Holm (1979)
Observations: 99

Result from sensitivity analysis correlation is close to the original correlation.